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Method of determining the degree 
of liquid aeration in a variable capacity 
displacement pump 


Zygmunt Paszota, Prof. 
Gdansk University of Technology, Poland 


ABSTRACT 


which the increase AMp,, App; Py pev 


indicated torque M, 
g Pi|App, = Py 4pgy 


The increase AM, 


Pm|App, = Py 4pgy 


The Author concludes, that there is a possibility of determining a concrete value of the 
liquid aeration coefficient e during the pump operation by finding such value of e with 
of torque of mechanical losses is proportional to the 
determined with a fixed value App, = cte of increase of 
pressure in the pump working chambers. The fixed value App; assumed in searching the liquid 
aeration coefficient e equals to the nominal pump operation pressure p, (App, = cte = p,). 
of torque of mechanical losses with a fixed value of App; 


(App, = cte) is proportional to the pump geometrical working capacity qp,,, therefore: only with taking 


into account the aeration coefficient e of liquid displaced by the pump the relation AM, 


m|dp Pi Pn IPgy ~ dpgy 


can be obtained from tests. The method, proposed by the Author, of determining the working liquid aeration 


coefficient e, is presented in this paper and has been practically applied for the first time by Jan Koralewski 


in his investigations of the influence of viscosity and compressibility of aerated hydraulic oil on volumetric 
and mechanical losses in a pump of HYDROMATIK A7V.58.1.R.P.F-00 type [8, 9]. 


Key words: hydrostatic drive; variable capacity displacement pump; liquid aeration; 
method of determining the degree of liquid aeration 


INTRODUCTION 


In references [1- 4] the Author evaluated the effect of 
working liquid compressibility on the picture of volumetric 
and mechanical losses in a high-pressure variable capacity 
displacement pump. The considerations are based on the 
assumptions made in the, developed by the Author, theoretical 
and mathematical models of torque of mechanical losses in 
a pump used in a hydrostatic drive [5 —7]. It is assumed in 
the models, that increase Mbisp,, a apay of torque of mechanical 
losses in the pump „working chambers - shaft” assembly is 
proportional to torque M,,; indicated in the pump working 
chambers: 


AM 


In references [1— 4] the Author introduced also a working 
liquid compressibility coefficient kiep, * Coefficient Kip, 
determines the decrease, as an effect of liquid compressibility, 
of active volume of working liquid displaced during one pump 
shaft revolution with the increase Apy, = p, of pressure in the 


Pm|App; ~ Mp; 


working chambers equal to the pump nominal pressure p, , 
compared with active volume equal to theoretical working 
capacity qp,or geometrical working capacity qp,, per one shaft 
revolution, determined with the increase Ap,, of pressure in the 
working chambers equal to zero — Apy, = 0: 


Art dptlAppi=Pa 
dpt 


k 


Iclp, = 


and 


dPgv7 dpgy|Appi= =Pn_ bp qpr bp dptlAppi=Pn 


pn bp dp 


A Pey 

The Author has also concluded that it is possible to evaluate 
the effect of liquid compressibility coefficient k,,,, on the value 
of increase AMp,, APpi = Pie Apay of torque of mechanical losses in 
the pump „working chambers - shaft” assembly and to evaluate 
the effect of kip coefficient on the value of the coefficient 
of volumetric losses in the pump working chambers due to 
leakage. 
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Searching the value of liquid compressibility coefficient 
kiep, , which with the increase App; of pressure in working 
chambers, equal to the pump nominal pressure p, , will cause 
an increase AM ES of torque of mechanical losses 
proportional to qp,, , i.e. to the indicated torque Mp, App; => Spey 
the Author determined, in the pump HYDROMATI 
ATV.58.1.R.P.F.00 type, tested by Jan Koralewski in his doctor 
dissertation [8] the approximate value of the oil compressibility 
coefficient during the tests as equal to kyga mpra = 0.030. 

Taking into account the working liquid compressibility 
evaluated by the coefficient k,,)35 mra = 0.030. Author determined 
approximate values of new coefficients of volumetric and 
mechanical losses in the tested pump. 

The Author concludes, that there is a possibility of 
determining a concrete value of the liquid aeration coefficient ¢ 
during the pump operation by finding such value of ¢ with which 
the increase AM,,, Apo: Pr don of torque of mechanical losses is 
proportional to the indicated torque Mp, OEI determined 
with a fixed value App; = cte of increase of pressure in the pump 
working chambers. The fixed value Ap,; assumed in searching 
the liquid aeration coefficient e equals to the nominal pump 
operation pressure p, (App, = cte = p,). 

The increase AM,,, APpi= Pi dp, of torque of mechanical losses 
with a fixed value of Apy, (Apy,= cte) is proportional to the pump 
geometrical working capacity qp,, , therefore: 


Pm|App; = P 


only with taking into account the aeration 
coefficient ¢ of liquid displaced by the pump the relation 


AMp mip, =P; apg ~ Ape CAN be obtained from tests. 


The method, proposed by the Author, of determining 
the working liquid aeration coefficient e, is presented in this 
paper and has been practically applied for the first time by Jan 
Koralewski in his investigations of the influence of viscosity 
and compressibility of aerated hydraulic oil on volumetric 
and mechanical losses in a pump of HYDROMATIK 
ATV.58.1.R.P.F.00 type [8, 9]. 


COMPRESSIBILITY OF WORKING LIQUID 
IN THE PUMP 


Compressibility of liquid at a given temperature is 
characterized by variation of its mass density p as a function 
of pressure p. In order to facilitate the calculations, the curve 
of variation p = f(p) is represented by an approximate algebraic 
relation. Most often a linear approximation is used: 


EE a) 
p B 

It may be said, that relation (1) defines a modulus B of the 
liquid volume elasticity valid for a certain temperature and 
for a certain pressure. 

Numerical indications regarding the value of B of hydraulic 
oils are the following [10]: 
— atthe normal temperature (20 °C), close to B = 1500 MPa 

for the used hydraulic oils, 
— B increases with the pressure (by about 1% every 2 MPa 

up to 20 MPa for the used oils (a, = 0.005/1 MPa)), 
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— Bdecreases when the temperature increases (about 1% every 
2 °C up to 100 °C for the used oils (as= - 0.005/1 °C)). 


In working chambers of the tested piston pump [8, 9], during 
their connection with the inlet channel, was slight overpressure 
Ppi; ~ 0.05 MPa (i.e. absolute pressure p,,,, ~ 0.15 MPa). Let's 
assume that the value of modulus B of the hydraulic oil volume 
elasticity, at the temperature 0 = 20 °C, equals to: 


Binetti wa = 1500 MPa (2) 


Therefore, the dependence of modulus B of oil on the 
increase App; of pressure in the working chambers and on 
the increase Ad of oil temperature may be described by the 
expression: 


B = B ppa = 0.15 Mpa; o=20°c (1 + a,App; + a940) (3) 


Modulus B of hydraulic oil volume elasticity decreases 
very quickly when oil is aerated, i.e. when the oil aeration 
coefficient e is greater than zero (£ > 0). 

The oil aeration coefficient ¢ is the ratio of volume V, of 
air to volume V, = V, + V, of mixture of oil volume V, and 
air volume V, (e = V,/V, = V,/(V, + V,)). The oil aeration 
coefficient e is determined at the absolute pressure py,,, in the 
pump working chambers during their connection with the inlet 
channel. 

Let's suppose, that volume V, of aerated oil in the pump 
working chambers, at initial absolute pressure py,, in the 
chambers (Fig. 6), contains a volume of air equal V, = V, and 
a volume of oil equal V, = (1 —s)V,. 

An increase App; of pressure in the pump working chambers 
causes a decrease of the oil and air mixture volume by the 
value AV (assuming a hypothesis of compression of air 
pV, = cte) equal to: 


V 
8 APpi (4) 
Pplia + Api 
Ifthe aeration coefficient e is small, which is a general case, 
V, is close to V, . Therefore it can be written: 


AV=AV,+AV,= = App; + 


AV=V, $ $ i] App, (5) 
B Pptia + APpi 

Therefore, with the oil aeration coefficient e greater 
than zero (e > 0), modulus B of oil volume elasticity must 
be replaced by modulus B’ of volume elasticity defined by 
the relation: 


1 l € 
SE e 
B B Pptia + APpi 
or, in the conditions of changing pressure and temperature of 
the aerated oil, by the relation: 


(6) 


1 1 
— = + 
B | 0.15 MPa, 9-20°c | + ApAP p; + 2509) 
€ (7) 
— 
Pplia + APpi 


B |ppyia=0-15 MPa,9=20°C = 1500 MPa 


a= 0.005/1 MPa aş = -0.005/1 °C 


0 5 10 15 
Indicated increase App, ofpressure [MPa] 


20 25 30 35 


Fig. 1. Modulus B of volume elasticity of non-aerated hydraulic oil (e = 0) and modulus B'of aerated oil (e > 0) as relations dependent on indicated increase 
App; of pressure in the pump working chambers, with limit values © = 20 °C (continuous line) and Bd = 68 °C (dashed line) of hydraulic oil temperature 
adopted during the investigations [8, 9]. It was assumed that modulus B of oil volume elasticity at absolute pressure pp, % 0.15 MPa in the pump working 
chambers during their connection with the inlet channel and at oil temperature 8 = 20 °C is equal to B = 1500 MPa. Also assumed was the value of 
coefficient a, = 0.005/1 MPa of the change of modulus B of oil due to increase App, of pressure in the working channels and coefficient a, = - 0.005/1 °C of 
the change of modulus B due to change of oil temperature 8 


Fig. 1 presents the modulus B of volume elasticity of 
non-aerated hydraulic oil (e = 0) and modulus B’ of aerated 
oil (e > 0) as relations dependent on indicated increase App, 
of pressure in the pump working chambers, with limit values 
0 = 20 °C and 0 = 68 °C of oil temperature adopted during the 
investigations [8, 9]. 

The displacement pump with variable geometric 
working capacity qp,, per one shaft revolution is tested with 
different fixed values of qp,, . 

Variable (set during the test) geometrical working capacity 
qpg, Of working chambers per one shaft revolution equals to 
the difference between maximum chambers capacity (capacity 
increasing to that maximum value during connection of the 
chamber with the pump inlet channel) and minimum chamber 
capacity (capacity decreasing to that value during connection 
of the chambers with the pump outlet (discharge) channel). The 
initial volume V, of oil (Fig. 6), which is compressed due to 
increase App, of pressure in the pump chambers, corresponding 
to setting qp,, of the variable geometrical working capacity, is 
in a variable capacity pump equal to: 


Vo = 0.5qp + 0.5qpgy (8) 


When the variable (set) geometrical working capacity qp,, 
reaches the maximum value equal to the pump theoretical 
working capacity qp; (qpgy = qr), volume V, of compressed oil 
is equal to: 


Vo = 0.5qp, +0.5qp, = qp: (9) 


The change AV of liquid volume (Fig. 6), due to the liquid 
compression under the increase App; of pressure in the pump 
chambers, is equal to the loss qp, of pump capacity per one 
shaft revolution: 


AV s dpve 


The loss qp,. of pump capacity per one shaft revolution 
(Fig. 2), resulting from the compression of non-aerated (or 
aerated) oil at the setting qp,, of geometrical variable working 
capacity, is determined (in reference to (5) and (6)) by the 
expression: 


(10) 


= (0.5qp, + 0.54p, )APp; 


APve B a 1) 
and with qp,, = dp; , by expression: 
qp: AP pi 12 
dpyc = M a ( ) 


and, after replacing 1/B’ by expression (7), by the formula: 
dpvc= (0.5qp, + 0.5qpgy y 


1 
+ 
(1+a,App; +asA9) (13) 


B 
|bpiia= 0.15 MPa,9=20°C 
í App; 


€ 


+ A 
Pplia + APpi 
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=1.000 


B = 1500 MPa 
a, = 0.005/MPa 


25 ay = -0.005/1°C 
i e=0.0135 


[cm?/rev] 


2.0 


Volumetric loss qp,, due to oil compressibility 


0.0 
0 5 10 15 


Coefficient bp of pump capacity Vn = 35mm’s! 


20 25 30 35 


Indicated increase App,ofpressure [MPa] 


Fig. 2. Loss qp,. of pump capacity during one pump shaft revolution due to compressibility of aerated (e = 0.0135) liquid, decreasing the active volume 
of liquid displaced by the pump compared with the theoretical working capacity qp, (bp = 1) or geometrical working capacity qpe, (0 < bp < 1) (pump of 
HYDROMATIK A7V.DR.1.R.P.F:00 type) [8, 9] 


and also, with qp,,= qr by the formula: 


dpvc "= Apt i 


1 


+ 
B Pp ja ¥0-15 MPa, 9=20°C (1+ a,App; +a5A9) 


— È _ 
Pplia + APp; 


Fig. 2 presents an example (with assumed oil aeration 
coefficient e = 0.0135) of the results of calculations of the 
loss qpy. = f(Apy.) of the tested pump capacity per one shaft 
revolution, taking into account formula (13) for setting qp,, of 
the geometrical variable working capacity and formula (14) for 
the maximum setting qp,, = qr: , 1.€. the case of pump theoretical 
working capacity. 

The change of qp, as a function of indicated increase App; 
of pressure in the working chambers, presented in Fig. 2, takes 
into account the influence of changing volumes V, (Fig. 6) of 
compressed liquid in the working chambers, the changes being 
an effect of the principle of operation of a variable capacity qpy, 
per one shaft revolution displacement pump (with variable b, 
coefficient). 

Loss qp,, of pump capacity during one shaft revolution due 
to the liquid compressibility decreases the active liquid volume 
displaced by the pump compared with the theoretical working 
capacity q» or geometrical variable working capacity qp,, 
(determined at Ap,, = 0). This fact should be taken into account 


+ 
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in the evaluation of intensity qp, = Qp,/np of volumetric losses in 
the working chambers as well as in the evaluation of increase 
AMbmjappi of the torque of mechanical losses in the „working 
chambers - shaft” assembly, the losses caused by the increase 
App; of pressure in the pump working chambers with determined 
values of the chamber geometrical working capacity qpy, . 

In the considerations, notions have been introduced 
of theoretical active working capacity qp, App; = Pa and of 
geometrical active working capacity qp,, Appi= p, 2S capacities 
the pump has in the working chambers at the increase App; 
of pressure in the chambers equal to the nominal pressure 
Pan of system operation. The active working capacities qpy Appi = 
p, and qpyy App; =p, CAN be determined from the equations: 


(15) 


dpilApp; =P = dp: m dpvejp, 


dPsv|App; =P, = dpgv dpvcjp, (16) 

Also a notion has been introduced of the coefficient k, 
of working liquid compressibility in the pump. 
Coefficient Kip, of the working liquid compressibility in 
the pump determines the degree of decrease, as an effect 
of liquid compressibility (without taking into account the 
chamber leakage), of the liquid active volume displaced by 
the pump during one shaft revolution at the increase App, 
of pressure in the pump working chambers equal to the 
nominal pressure p, of system operation, compared with the 
active volume displaced by the pump at App; = 0. Coefficient 
kip, is defined by the formulae: 


lelp,, 


qre 7 dPt|Appi=PA dPve|p, 
ie|p, = = (17) 
qrt qPt 
or: 
e dpgv m dpgv|App, =Pa _ dpvejp, (18) 
telpa q = q 
Pgv Pvg 


The knowledge of coefficient k, of the liquid 
compressibility in the pump allows to evaluate numerically 
the subdivision of volumetric losses in the pump into the 
losses due to oil leakage in the working chambers and losses 
due to liquid compressibility. 

In a variable capacity pump operating with setting qp,, of 
geometrical variable working capacity (determined at Apy, = 0), 
the coefficient k,,, is described (in reference to (13) and (18)) 
by the formula: 

B A PvdAppi=pn 0.5qp, + 0.54 poy 
lepn = i 


Iclp, 


dpgv dPgv 
| + 
(19) 
B ia s015 MPa,9=20 °C (+ ApPn + ag AS) 
: Pn 
E€ 
+ — 
Pplia + Pn 


and with qp,, = qr; (in reference to (14) and (17)) by the 
formula: 


_ dPvdApp;=P, _ 


lep dpi 


k 


1 


+ 
(1+a,Pn +a, A9) > 


B 
[priia~0.15 MPa,9=20°C 
= Pn 


e 


PPlia + Pn 

Therefore, in a displacement pump, operating at theoretical 
working capacity qp per one shaft revolution, the working 
liquid compressibility coefficient ki (formula (20)) is 
a combined effect of modulus B of hydraulic oil volume 
elasticity, oil aeration coefficient e and also liquid temperature 
Ò (increase AS to the reference temperature 0 = 20 °C) and 
absolute pressure Ppi in the working chambers during their 
connection with the inlet channel, as well as the system 
nominal pressure p, . 

In the same displacement pump operating with variable 
capacity qp,, per one shaft revolution, the value of working 
liquid compressibility coefficient Kip (formula (19)) 
increases in comparison with the value Ex during the pump 
operation with the theoretical working capacity qp, . This 
is an effect of the increased ratio of initial liquid volume 
(V, in Fig. 6) subjected to compression, i.e. the volume 
(0.5qp, + 0.5qp,,) (formula (8)), to the set working volume 
dps, - Therefore, decreasing the qy,, setting causes in 
a variable capacity displacement pump an increase of 
Kip coefficient (formula (19)). 


+ 


Ielp, 


IMPORTANCE OF THE ACCURACY OF 
EVALUATION OF qp; AND q,,, FOR THE 
ACCURACY OF EVALUATION OF THE 
INTENSITY OF VOLUMETRIC LOSSES 
AND TORQUE OF MECHANICAL LOSSES 
IN THE PUMP 


It is important, particularly during investigation of a pump 
of variable capacity per one shaft revolution, to determine 
precisely the theoretical capacity qp, per one pump shaft 
revolution and geometrical capacities qp,, per one pump 
shaft revolution. Geometrical capacities qp,, change in 
the 0 < qp,, < qp, range and the corresponding coefficients 
bp = qp.//qp. of the pump capacity change in the 0 <b, <1 
range. Therefore, precise evaluation of the coefficient 
bp = qpgv/dp. depends on the precision of evaluation of qpy, 
and qp, . 

The pump theoretical working capacity qp, and geometrical 
working capacities qp,, are determined at the indicated increase 
App; of pressure in the pump working chambers equal to 
zero (App; = 0); the working capacities are determined by 
approximation at the point App; = 0 of the line qp = Q,/n, = 
= f(Ap,,) describing, at a pump constant setting (but the value 
of b, coefficient unknown exactly), the capacity qp displaced 
during one pump shaft revolution as dependent on the 
indicated increase App; . The line qp = f(Ap,,) runs through the 
measurement points obtained in the investigation. 

Fig. 3 presents examples of dependence qp = f(Apy,) of 
capacity qp per one shaft revolution of an axial piston pump 
on the indicated increase App; of pressure in the working 
chambers at the coefficient b, = 1 of pump capacity per one shaft 
revolution [8, 9]. These are examples of searching the theoretical 
working capacity qp, per one pump shaft revolution and also the 
evaluation of subdivision of the intensity qp, of volumetric losses 
per one shaft revolution into the volumetric loss qp,„ due to oil 
leakage in the working chambers and volumetric loss qp,, due 
to compressibility of non-aerated (or aerated) oil. 

Loss qp,. = f(App,) per one shaft revolution, determined from 
formula (13), due to liquid compressibility, occurring at setting 
pg, Of the pump variable geometrical working capacity (or 
from formula (14) at setting q,, of the pump theoretical working 
capacity) is added to capacity qp= f(Ap,,) per one shaft revolution 
shown by the line running through the test measurement points. 
The result of adding qp,. = f(App;) to qp = f(Apy,) is pump capacity 
dp without compressibility ~ {(APp;) as a difference between qp,, (Or qp,) 
and the volumetric loss qp,, due to oil leakage (independent of 
the liquid compressibility): 

(Gp without compressibility = Arve T dp) = (Apr) (21) 
(Gp without compressibility APgv ( or dp.) — Gp) = f(App;) (22) 
Approximation of the qp without GE = f(Ap;, line at 


App, = 0 allows to determine the qp,, (Or p,): 


dp without compressibility |ApPi = 0 = dpgv (or dpi) (23) 


As shown in Fig. 3, the pump theoretical working capacities 
qp. , at point App; = 0 of the qp = f(Ap,,) line, obtained from tests 
and including also the liquid compressibility, as well as of the 
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bp= 1.000 v/v, = 1.000 


qr: at £= 0.016 


Volumetric loss qp,, dueto 
compressibility of non- aerated 
(e= 0) or aerated (£> 0) oil 


Capacity qpper one shaft revolution [cm?/rev] 
w 
> 
in 


w 
> 
o 


w 
> 
in 


53.0 


B |ppqia 70-15 MPa,9=20°C 7 1500 MPa ap= 0.005/1 MPa ag = -0.005/1 °C 


Capacity qp per one pump shaft i 
revolution marked by measurement 
points obtained from tests 


Va = 35 mms! 


qr at £ = 0.000 


Volumetric loss qp,; 
due to oil leakage 


/ ł 
Capacity qp without compressibility = | 
taken into account 
m www www | mm] 
20 25 30 35 


Indicated increase App; ofpressure [MPa] 


Fig. 3. Dependence of pump capacity qp per one shaft revolution on the indicated increase App, of pressure in the working chambers, at the coefficient 
bp = 1 of pump capacity; the val-ues qp,, of geomet rical working volume and qPt of theoretical working volume per one shaft revolution (determined 
at App; = 0) and subdivision of the intensity qp, = py + qp, of volu-metric losses per one shaft revolution into volumetric loss qp,, due to oil leakage in the 
cham-bers and volumetric loss qp,. due to compressibility of non-aerated (or aerated) oil dependent on the value of oil aeration coefficient € 
(e = 0 to 0.016); viscosity coefficient v/vn = 1, oil temperature © = 43 °C (pump of the HYDROMATIK A7V.DR.1.R.P_F.00 type) [8, 9] 


(4p without compressibility = Gpve + dp) = f(APp)) line taking into account 
compressibility of non-aerated (e = 0) oil, have practically 
the same value qp, = 58.9 cm3/rev. Approximation of the 
(4p without compressibility = pve + dr) = f(App:) line at App, = O point, 
taking into account the aerated oil compressibility, shows an 
increase of qp, value practically proportional to oil aeration 
coefficient e. This is clearly presented in Fig. 4. For example, 
theoretical working capacity with e = 0.0135, reaches the value 
dp, = 59.57 cm3/rev. 

Fig. 5 presents the subdivision of volumetric losses qp, = 
= f(Ap,.) in the tested pump into loss qpe = f(Ap,;) due to 
liquid compressibility and loss qp, = f(App;) due to oil leakage 
at different values of liquid aeration coefficient e, with the 
theoretical working capacity q,, per one pump shaft revolution. 
The lines of qp, = f(Ap,;) loss due to oil leakage do not change at 
different oil aeration coefficient e value, but lines qp,, = f(App.) 
due to oil compressibility differ clearly, as well as lines (qp, = 
= dpu + dp.) = f(APp.) of volumetric losses qp, = f(Ap,;) as a sum 
of qp, = f(Ap,;) of loss due to leakage and qy,, = f(Ap,;) of loss 
due to liquid compressibility. 

During a careful investigation of a displacement pump 
pressing non-aerated working liquid, the accuracy of 
determining the theoretical working volume qp, and geometrical 
working volumes qp,, is of an order of one per mille. Accuracy 
of estimation of the value of coefficients b, of pump capacity 
is then also very high. 

The accuracy of evaluation of qp, and qp,, is significantly 
worse when the working liquid is aerated. This is an effect of 
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bp = 1.000 Va = 35 mms"! 
B |ppq.q=0:15 MPa,8=20°C = 1500MPa a„=0.005/1 MPa ag = 0.005/1 °C 
60.0 


59.5 
y= 55.451x + 58.907 


[cm*/rev] 


59.0 


Theoretical working capacity qp, per one shaft revolution 


58.5 
0.000 


0.016 


0.004 0.008 0.012 0.020 


Oil aeration coefficient e 


Fig. 4. Effect of evaluation of theoretical working capacity qp, per one pump 
shaft revolution resulting from assumption of aeration coefficient e of the 
pump displaced oil; evaluation qp, (Fig. 3) is a result of approximation, at 
App; = 0. of the relation of pump capacity qp per one shaft revolution to the 
indicated increase App; of pressure in the working chambers taking into 
account the aerated oil compressibility (at a given oil aeration coefficient e) 
(pump HYDROMATIK A7V.DR.1 .R.P. F00 type) [8, 9] 


the fact, that the aerated liquid in the working chambers filled 
during their connection with the low-pressure inlet channel 
decreases its volume because of great compressibility of non- 


Va = 35 mms! 


bp = 1.000 v/v, = 1.000 


Oil aeration coefficient e 


e=0.000 €= 0.004 e=0.008 £:=0.012 2¢=0.0135 €= 0.016 
7.0 
6.5 B = 1500 MPa 
: a, =0.005/1 MPa 


ag = -0.005/1 °C 


6.0 Volumetric loss qpy = qpv + Gpve 
5.5 
Volumetric loss qpy 
5.0 due to oil leakage , 
4.5 \ 
4.0 Volumetric loss qpy, due to \ 


compressibility of non-aerated \ 
(e = 0) or aerated (e > 0) oil 


we 
w 


Volumetric loss qp,, Gpycs qpu per one shaft revolution [cm*/rev] 
w 
o 


35 


Indicated increase App, ofpressure [MPa] 


Fig. 5. Subdivision of volumetric loss qp, = J(App,) in the pump into loss qp, 
= f(App,) due to oil compressibility and loss dp, = f(App;) due to oil leakage 
at different values of oil aeration coefficient e and value v/v, = 1 of oil 
viscosity coefficient in the tested pump with theoretical working capacity qp, 
(bp = 1) (pump HYDROMATIK A7V.DR.1.R.P.F.00 type) [8, 9] 


dissolved air in the liquid, after connection of the working 
chambers with the discharge channel where pressure may be 
only a little higher than in the inlet channel. 

Without knowledge of the coefficient e of aeration 
of the liquid flowing into the pump working chambers, 
it is impossible to determine the quantities qp, and qp,, 
precisely. 

At the same time, precise knowledge of qp, and qp,, is 
important in evaluation of the volumetric and mechanical 
losses in the pump. 

The intensity qp, = Q,,/n, of volumetric losses Qp, in the 
pump working chambers per one pump shaft revolution is 
evaluated as a difference between qp, (or qp,,) and qp determined 
during the investigation at changing values of the indicated 
increase Ap,, of pressure in the chambers. 

Increase AM;,, App; of the torque of mechanical losses in 
the pump „working chambers - shaft” assembly, compared 
with torque Mp, i of mechanical losses in the no-load pump, 
is an effect of increased friction forces in the assembly resulting 
from the influence of the torque M,, indicated in the working 
chambers upon the assembly and is proportional to M,, . 

Increase AMpmiapy, of the torque of mechanical losses in 
the pump ,,working chambers - shaft” assembly is determined 
during the investigations as a difference AM bmnp; = M, + 
— Mom AR between torque Mpm of losses in the assembly and 
torque Mp, APpi -o of losses in the assembly of a no-load pump. 


Torque M,, of losses is determined as a difference M,,, = 
= M, — M; between torque M, measured directly on the shaft 
and torque M,, indicated in the working chambers. Therefore, 
extremely important, for determination of torque M,,, of 
mechanical losses and of increase AM miap,, of torque of 
mechanical losses, is precision of determination of torque 
M,, indicated in the working chambers (defined by formulae 
(35) and (36)). 


WORK OF THE DELIVERY OF 
COMPRESSED WORKING LIQUID 
DURING ONE PUMP SHAFT REVOLUTION 
AND TORQUE INDICATED IN THE PUMP 
WORKING CHAMBERS 


In order to deliver compressed working liquid during one 
pump shaft revolution, a certain work E is required, which 
is a sum of: 

— the work of compression itself — E,, 
— the work of displacement at constant pressure — E, . 


Let’s calculate the theoretical values (efficiencies equal to 
1) of the two works. For that purpose, the pump delivering the 
compressed liquid is presented in a simplified way by a piston 
of cross-section S moving in a cylinder, which, by two valves 
R, and R, (acting as a distributer), can communicate with, 
respectively: 
— space filled with liquid at constant absolute pressure Ppiia , 
— volume C, filled with liquid at constant absolute pressure 


Pra. (Fig. 6). 


Piston in Fig. 6 performs stroke from position x, to the 
bottom of cylinder i.e. to position 0 


Do U 


0 X; Xo x 


Fig. 6. Work of the delivery of compressed working liquid during one 
displacement pump shaft revolution (simplified diagram) 


Initial position: piston at x, . Volume open to liquid at 
pressure pp, : Vo . R, opened, R, closed. 
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1st phase: R, is closed, piston pushed from x, to x, (volume V,). This is the point where liquid, closed in the cylinder, achieves 
the pressure pp>,, . The work delivered by the piston is the work of compression: 


X, M 
E=- fe — Ppiia ) Sdx =— [ (P= Pria) AV FM 
Xo Vo 


Work E, is represented by area ABFA. 
2nd phase: R, is opened and the piston is pushed (in the case of bp = qp,,/qp, = 1) from x, to 0. i.e., to the bottom of cylinder. 
The liquid is discharged into C, . The work delivered by the piston is the work of displacement: 


0 0 
E, =- [poi — Pplia )Sdx =- [Orza —Pptia AV = ( Ppzia™ PPlia) Vi (25) 


xy Vi 


Work E, is represented by area BCDFB. 

3rd phase: R, is closed, R, opened and we return to the initial position. This operation is carried out without work performed 
by the pump. 

The total work E = E, + E, is represented by the dashed areas in Fig. 6. 

One of the definitions of the modulus B of liquid volume elasticity is the following: 


AV Ap V 
— =-— where V =-—d 2 
y B B P (26) 
Therefore, the compression work is: 
Vi Po 2 
V V (p —p,) P2 
E, == —p,)dV= —p)—dp=—*—+ 
z fe pi) fe Doo h (27) 
0 


During operation of pump with theoretical capacity qp, per one shaft revolution (b, = 1), the change of V (Fig. 6) during 
compression work in relation to Vis small. The compression curve may be replaced by linear approximation and the V quantity 


in equation (27) by mean value Vn = (V, + V,)/2: 
E -Vn (2 =P (28) 
= 
B 2 
Therefore: 
E, = V; (Praia — Ppria) (ace. to (25)) 
and: 
V — AV 
E =E, +E, = (p, - pV, + 2a =P) =(P = P) Y + 
2B 2 
E = (P> m POV 
E = V m (Przia — Priia) (29) 
Formula (29) describing the work E may be replaced by the expression: 
AV AV 
E =v- (P P2ia™ Ppiia) -[v- = Appi (30) 


In a real displacement pump with variable capacity per one shaft revolution, with the geometrical working capacity setting 
pay » Work E performed by the pump in the working chambers during one shaft revolution (after replacing in formula (30) the 
original volume V, by volume qp,, , the change AV of liquid volume due to liquid compressibility by loss qp,, of pump capacity 
during one shaft revolution (formula (8) and the loss qp,. by formula (11)), is described by the expressions: 


q VC 
E= (ars = > App (31) 
and: 
E=41 : l + ——— [Apri > dpo, Ap (32) 
=0 ® Pi Pgv © F Pi 
ele Ppqi. 70.15 MPa,9 =20 (l T a APpi tag A9) Pplia + APpi i i | 
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and with qp,, = dp; (in reference to (12)) by the expressions: 


E= fan- EB App (33) 
and: 
1 1 E 
E=11-- +————- Apy AD»; 
2|B (1+a,App; + ag AS) Drin + App; Ppi¢ Apt OP pi (34) 


Q 
[Ppa 0.15 MPa,9 =20°C 


It has to be mentioned, that in formula (32), describing the work E performed by the pump with setting qp,, per one shaft 
revolution, the qp,, value is determined by approximation of the (qp + qp,.) = f(App,) line at point App, = 0. At the same time, the 
formula (32) contains the value q,,,. = f(Ap;,) describing a loss of capacity per one shaft revolution due to the liquid compressibility 
taking into account the change AV (Fig. 6) of liquid volume resulting from the variable capacity pump operation mode, i.e. the 
compressed volume V, (Fig. 6) equal to V, = 0.5qp, + 0.5qpz, . 

Torque M,,; indicated in the pump working chambers, corresponding to work E in the chambers during one shaft revolution, 
is then, with qp, setting, described by the formula: 


AD». 
Na = 1 PREZ. AAAA RSE: AAA App; dpgy OP pi (35) 
211 2 Bp, =0.15 MPa,9=20°C (l +a,APp; + ag A9)  Ppiia + APpi 211 
and with qpy, = qp. by the formula: 
AD». 
meh o OC 
211 2 Bp, =0.15 mPa,9=20°C (l +a App; + Ag A9)  Ppia + APpi 211 


METHOD OF DETERMINING THE WORKING LIQUID AERATION COEFFICIENT s 


As it so far has not been possible to determine the coefficient ¢ of aeration of the working liquid flowing into the pump, 
and therefore not possible to take the liquid compressibility into account, both at the small increase App; of pressure in the pump 
working chambers and in the full range of increase Ap,; of pressure — up to the hydrostatic drive system nominal pressure p, , 
the pictures of volumetric losses and mechanical losses in the pump, determined by the applied methods, are deformed. 
For instance, if the compressibility of working liquid, the liquid in fact aerated, is not taken into account, a picture of negative 
increase AMpmiapp; of torque of mechanical losses in the pump ,,working chambers - shaft” assembly is obtained as an effect of 
increased torque Mp, indicated in the working chambers due to the increase of geometrical working volume qp,, (of bp coefficient), 
which is an illogical result (Fig. 7). 

A method of determining the coefficient ¢ of working liquid aeration may be searching for value of ¢ which was used 
for determining the values qp,, of geometrical working capacities causing an increase Momispp; of torque of mechanical 
losses in the pump ,,working chambers - shaft” assembly proportional to torque M,, indicated in the working chambers 
(formula (35)), a torque Mp, resulting from q,,, and from ze at the constant value App, of indicated increase of pressure 
in the chambers. It is assumed, that during searching for qp,, and e, the increase Mp App; of the torque of mechanical losses is 
determined at a constant value of indicated increase App; of pressure in the pump working chambers equal to the system nominal 
pressure p, (App; = P,)- 

Therefore, it is assumed, that with a fixed value Apy, = p, of the indicated increase of pressure in the pump working 
chambers, the increasing torque M,, in the chambers (formula (35)), described by the expression: 
MAI t p E p ae Po aż 

2 s. <0.15 MPa 9 -20°C +a pPn + Ag A9)  Ppra + Pn 2II 
must be accompanied by proportional increase AM 
- shaft” assembly: 


PmiAPp;=P, of torque of mechanical losses in the pump ,,working chambers 


AMb mapp; = Pp? dpgy ~ Mbijappi = Pp dpgy G 8) 
or 
AMpmjappi= Pp! dpgv ~ dpgv (bp) (39) 
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With fixed values of B, a, , ay , Ù, pp,,, and p, , expressions 
(38) and (39) can be obtained only with one value £ of the 
aeration coefficient, which was assumed for determining of 
pe, Values and of the pump capacity bp coefficient values. 

Fig. 7 presents results of searching for the oil aeration 
coefficient e during pump tests [8, 9]. 

With the assumption of non-compressible (B = œ) and 
non-aerated (e = 0) liquid, i.e. with the assumption of the liquid 
compressibility coefficient Kup, = 0 the picture of AMp,, Apit = 
= f(Mpia ) relation has the form of a descending straight 
line: 

from value AMp,, = 1.87 Nm at qp,, = 0 

to value AMp,, = 0.53 Nm at qp,, = qp, - 


Ppi Pn 


With the assumption of compressible and non-aerated 
(e = 0) liquid, the picture of AM,,, Ant (Mp, Appi =p) relation 
has the form of an ascending straight line: 

from value AMp,, = 1.86 Nm at qp,, = 0 

to value AMp,, = 2.79 Nm at qp,, = qr: - 


With the assumption of compressible and aerated (z = 0.008) 
liquid, the picture of AM,,,, bp; “Py f(Mp, Appi -p relation has 
the form of an ascending straight line: 

from value AMp,, = 0.76 Nm at qp,, = 0 


to value  AMp„= 1.77 Nm at qp,, = qp: - 


With the assumption of compressible and aerated (z = 0.016) 
liquid, the picture of AMpmiap,, = p, = f(Mpiap.,. = p) relation has 
the form of an ascending straight line: 

from value AM,, = — 0.35 Nm at qp,, = 0 


to value AMpm = 0.74 Nm at dp, = dp - 
B=1500 MPa 
KOLE EEG oc eee 
a 
g | 
2 | 
B3 a |--A AC O OI 
£ | e=0.000 + 
=) y = 0.9327x + 1.8617 | + 
v 2 EE Ei iacansaafencsncsnsseeanecssensancenennscnscscsncnnssnsssesescafesesncnnecenssnsnnees. 
S 
5 + B = œ MPa 
fra y = -1.3385x + 1.8663 + 
© 
e=0.008 ° 

È 

SS 

£ 


y = 1.0106x + 0.7632 


y= 1.0259x 


y = 1.0899x - 0.3512 


Increase AMpm| App; 


0.0 0.2 0.4 


ap = 0.005/1 MPa 


With the oil aeration coefficient e = 0.0135, the picture 
AMP Appi =P dpgy f(Mp, AE apg) has the form ofan ascending 
straight line: 

from value AMp,, = 0 at qpy, = 0 

to value AMp,, = 1.03 Nm at qp,, = qp: - 


Figure 8 presents, in reference to results presented in 
figure 7, a linear relation of the increase AMomiapp; =p; apg =° of 
torque of mechanical losses to the assumed value of the oil 
aeration coefficient e. The relation shown in figure 8 allows to 
find with high accuracy the value of oil aeration coefficient 
e, with which the increase AM pnjap; = Pa dpgy =O of the torque 
of mechanical losses at qp,, = 0 (bp = 0) equals zero: 


Pm|Appj = Py} dpgv "3 (40) 

The oil aeration coefficient ¢ during the pump 
testing (pump HYDROMATIK A7V.DR.1.R.P.F.00 type) 
corresponding to the situation described by formula (40), 


had the value ¢ = 0.0135 [8, 9]. 


CONCLUSIONS 


1. The Author concludes, that there is a possibility of 
determining a concrete value of the liquid aeration 
coefficient e during the pump operation by finding such 
value ofe with which the increase AM,,, Api Py apay of torque 
of mechanical losses is proportional to the indicated torque 
AMp App; = Pa Apg determined with a fixed value App; = cte of 

increase of pressure in the pump working chambers. 


ag = -0.005/1 °C 


hw ww w w a wwa ść 


0.6 0.8 1.0 


Pump capacity coefficient bp 


Fig. 7. Picture of the relations of increase AM,,,, Pp; = 32 MPa; dp, 


of torque of mechanical losses in the pump „working chambers - shaft” assembly (pump 


HYDROMATIK A7V.DR.1.R.P.F.00 type) to the geometrical working capacity qp,, (bp coefficient) with assumed values of modulus B of hydraulic oil elasticity 
and oil aeration coefficient e; the line AMpmy, - 32 urarg, Corresponding to e = 0.0135 is an effect of the straight line picture presented in Fig. 8 [8, 9] 
Pi” "Igy 
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—0.020 


Qil aeration coefficient 


05 “00 0.5 


Increase AMpm | App;= Pr. qpgv 


B=1500 MPa a, = 0.005/1 MPa 


ag = -0.005/1 °C 


| y=-0.0072x+0.0135 


1.0 1.5 2.0 


=o Of torque of mechanical losses [Nm] 


Fig. 8. Straight line relation of the oil aeration coefficient e to the assumed increase AM,,, dp pj =32 MPa: dp. = 0 of torque of mechanical losses 
(pump HYDROMATIK A7V.DR.1.R.P.F'00 type) [8, 9] 


2. The fixed value App; assumed in searching the liquid aeration 


coefficient e equals to the nominal pump operation pressure 
Pn (App; = cte = p,). 
The increase AM,,, AbD Sia of torque of mechanical losses 
with a fixed value of App; (App, = cte) is proportional to the 
pump geometrical working capacity qp,, , therefore: 
only with taking into account the aeration coefficient € 
of liquid displaced by the pump the relation 

AMp pp = Pyidpgy 7 Pe CAN be obtained from tests. 
The method, proposed by the Author, of determining the 
working liquid aeration coefficient e, is presented in this 
paper and has been practically applied for the first time 
by Jan Koralewski in his investigations of the influence 
of viscosity and compressibility of aerated hydraulic 
oil on volumetric and mechanical losses in a pump of 
HYDROMATIK A7V.58.1.R.P.F.00 type [8, 9]. 
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with volumetric control 
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ABSTRACT 


p 


There are uninvestigated areas connected with behaviour ofelements in hydraulic systems 
with different structures. Unawareness of proportions of the energy, volumetric, pressure 
and mechanical losses in elements is often the case. Problems connected with energy 
efficiency are essential for improvement of functionality and quality of hydrostatic drive 
systems, characterised by unquestioned advantages but also by relatively low efficiency 
in comparison with other types of drive. Energy efficiency of hydrostatic transmissions, 
a particularly those with throttling control of the motor speed, and also efficiency of the 


hydraulic servo-mechanism systems may be in fact higher than the values most often quoted in publications 

on the subject. Possibility of calculating the real value of the hydraulic system overall efficiency as 

a function of many parameters influencing it, becomes a tool of complete evaluation of the designed system 

quality. The paper compares efficiencies of systems with cylinder proportional control and efficiency of the 

system volumetric control by a variable capacity pump. Presented are also two schematic diagrams of the 

investigated hydrostatic systems, their principle of operation and problems of studying losses in elements 
and energy efficiency of systems consisting of a feed assembly, control set and cylinder. 


Key words: hydrostatic system; control structures; proportional 
directional valve; cylinder; energy losses; energy efficiency 


INTRODUCTION 


In searching for the energy saving solutions, computer-aided 
methods of calculating the energy efficiency of systems have 
been developed and improved. 

The hydrostatic systems play a very important role in 
modern machines. Great number of the nowaday constructed 
machines have more or less developed hydrostatic or electric- 
hydrostatic drive systems and in many cases those systems are 
the most important parts of the machines. Component elements 
— hydraulic linear motors (cylinders) — are widely applied 
in machines used on land and aboard ships. Unquestioned 
advantages of cylinders are: capability of performing the 
translational motion, reliability, simple construction, the 
effective force to weight ratio. 

The required speed v,, and load F,, of the driven machine 
are a result of its operation cycle and tasks to be performed. The 
driven machine current speed and load values are independent 
of the type and structure of the machine driving system. 

The current speed and load of the hydrostatic system driven 
machine have a direct or indirect impact on the mechanical, 
volumetric and pressure losses in the hydraulic motor, pump 
and other elements of the system with a given motor speed 
control structure, the losses resulting also from the hydraulic 
oil viscosity. 
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If, as an effect of the increasing, required by the driven 
machine, hydraulic motor operating speed v,, or increasing, 
required by the machine, motor load F,, , and also as an effect 
of the mechanical, volumetric and pressure losses in the 
hydrostatic driving system elements, the maximum driving 
system capability, determined by the maximum pump capacity 
Qpmax Or MAXIMUM pressure Ppzmax 1 the pump discharge conduit 
limited to the system nominal pressure p, , is used up, then 
further increase of vy, or Fy will not be possible [4--9]. 

The maximum pump capacity Q,,,,, is less than its 
theoretical capacity Qp. The pump theoretical capacity Qp, is 
a product of the theoretical capacity qp, per one pump shaft 
revolution and the no-load pump shaft speed n,,. The pump 
Qpmax Capacity results from the loaded pump speed np, lower 
than the np speed. At the same time, volumetric losses occur 
in the pump [4-9]. 

The system nominal pressure p, is a maximum permissible 
pressure pp, for its continuous operation, determined in the 
pump discharge conduit. 

Maximum values of the hydrostatic drive system hydraulic 
motor speed Vyuma and load Fuma are limited by the maximum 
pump capacity Q,,,,, and the system (pump) nominal pressure 
p, , and also by the mechanical, volumetric and pressure 
losses in other system elements, which are also an effect of 
the working liquid viscosity [4+9]. 


Fig. 1. Diagram of the investigated system fed at constant 
pressure — a p = cte structure [1+3] 


The most often used proportional control system of a linear 
hydraulic motor is a system (Fig. 1), where the proportional 
directional valve is fed by a constant capacity pump cooperating 
with an overflow valve stabilizing a constant feed pressure level 
(p = cte). Such system achieves high energy efficiency, close 
to the efficiency of a system without throttling control, only 
at a point of maximum values of the motor load coefficient My 
and speed coefficient wy. With decreasing motor load, and 
particularly with simultaneous decreasing motor speed, the 
system efficiency 7 decreases rapidly [1+3]. 


Apwi=Pui ~ Puzi 


a 


Fum 
Fui = Fm +Fum / 


Fig. 2. Diagram of the investigated system with proportional directional 
valve fed by a constant capacity pump cooperating with an overflow valve 
controlled in a variable pressure system — p = var [1+3] 


There are possibilities of reducing the energy losses in 
elements (pump, throttling control assembly and hydraulic 
motor, particularly the linear motor) of a system with 
proportional control, i.e. possibilities of increasing the energy 
efficiency of a system with throttling valve. 

The hydrostatic system of drive and linear hydraulic motor 
proportional control may be fed by a constant capacity pump 
cooperating with an overflow valve stabilizing the proportional 
directional valve feed pressure at the nominal pressure level 
(Fig. 1), or by a pump cooperating with an overflow valve 
controlled by the receiver inlet pressure. The variable pressure 
(p = var) (Fig. 2) system makes it possible to reduce losses in 
the pump, in the control assembly and in the linear hydraulic 
motor [1+3]. 

In a variable pressure (p = var) system, the structural 
pressure and volumetric losses in the throttling control 
assembly, mechanical losses in the cylinder and pump as well as 
volumetric losses in the pump can be significantly reduced. The 
mathematical description of losses and efficiency is presented 
in reference [1]. 

The hydrostatic system structure has an essential influence 
on the system efficiency. Its impact is most often considered 
with the assumption of an ideal pump and motor and also with 
supposition that energy losses in the real pump and motor 
will cause further proportional decrease of the system overall 
efficiency. However, the picture of mutual impact of losses 
in all the hydrostatic system elements appears much more 
complex [8]. 


IMPORTANCE OF THE INVESTIGATIONS 


Problems connected with energy efficiency are of basic 
importance for improvement of functionality and quality of 
the hydrostatic drive systems, characterised by unquestioned 
advantages but also by relatively low efficiency in comparison 
with other types of drive. Publications describing the 
influence of particular design and operating parameters on 
the hydrostatic system efficiency are valuable. They make it 
possible to work out system configurations with losses reduced 
to a minimum. 

Energy efficiency of hydrostatic transmissions, particularly 
those with throttling control of the motor speed, and also 
efficiency of the hydraulic servo-mechanism systems may be 
in fact higher than the values most often quoted in publications 
on the subject. Possibility of calculating the real value of 
the hydraulic system overall efficiency as a function of 
many parameters influencing it, becomes a tool of complete 
evaluation of the designed system quality. The capability of 
making such evaluation is important because the hydrostatic 
control systems are used in various machines and equipment, 
and also due to increasing power of the hydrostatic drive at the 
time of constantly increasing costs of energy generation [8]. 

In a system with too low efficiency, the load, first of all of 
the pump, increases, which leads to increased risk of failure 
and the necessary repair or replacement, as well as to a shorter 
service life. The too low system efficiency, most often resulting 
from intensive throttling of the stream of liquid, is also a source 
of rapid deterioration of operational characteristics, particularly 
the hydraulic oil lubricating properties, which is an effect, 
among other reasons, of too high temperature of the working 
liquid — the hydrostatic transmission power medium. 

Comparison of power losses in the elements provides 
information facilitating design of a new system. 

Comparison of energy balances from the point of view of 
the power of losses in different systems gives a broader view 
for choosing an optimum solution. 
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COMPARISON OF THE ENERGY 
EFFICIENCY OF VARIOUS SYSTEM 
VERSIONS 


Investigations of the efficiency of elements and systems, 
taking into account detail analysis of the sources of particular 
energy losses, may be included in the basic scope of research 
into the hydrostatic drive and control systems. 

The energy efficiency, one of the most important system 
characteristics, is defined as the ratio of, by the driven device 
currently demanded, useful power Py, of the hydraulic motor 
to the power Pp. , corresponding to Pyu , obtained by the pump 
on its shaft from the (electric or combustion) motor driving 
the pump. In case of improper choice of the system type, the 
consequence may be increased hydraulic oil temperature, 
i.e. decreased oil viscosity and, in turn, lower efficiency of 
the system elements, and also an impact on the system run 
characteristics. Therefore, the energy efficiency may be 
a decisive factor for usability of a system in a specific case. But 
its detailed analysis quite often leads to design improvements 
of the system elements. However, improving the quality 


of hydrostatic systems does not consist exclusively in the 
improvements of their elements [1]. 

Figures 3 and 4 present the overall efficiency n ofa constant 
pressure (p = cte) and a variable pressure (p = var) system with 
proportional control and a system with volumetric control 
by a variable capacity pump (Qp = var) as a function of the 
load coefficient My, at different values of the cylinder speed 
coefficient Oy. 

In the case of a system with volumetric control by 
a variable capacity pump (Q, = var), increasing the cylinder 
load coefficient My, causes rapid increase of the system overall 
efficiency n (Fig. 3). However, efficiency of structures with 
the series throttling control fed by a constant capacity pump 
is, with small value of the dy coefficient, distinctly lower than 
the volumetric control efficiency with the same value of Oy, 
because the structural losses are high. 

Increase of the cylinder speed causes a proportional 
increase of efficiency of the p = cte and p = var systems, but 
with an increase of the cylinder speed v,, the relative increase 
of efficiency of the system fed by a variable capacity pump is 
smaller (Fig. 3). 


The overall efficiency n 
of a constant pressure (p=cte) and a variable pressure (p=var) system and a system with volumetric control (Q.=var) 


- determined by simulation from experimentally obtained k; coefficients 
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Fig. 3. Relation of the overall efficiency n of a constant pressure (p = cte) and a variable pressure (p = var) system with proportional control and a system 


with volumetric control by a variable capacity pump (Q, = var) to the 


load coefficient My, at different values of the cylinder speed coefficient Oy 
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(efficiency determined by simulation from experimentally obtained k, coefficients ; the vy = 0.350m/s (Oy = 0.875) 
speed was the highest cylinder speed used during the tests) [1] 
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The overall efficiency n 
of a constant pressure (p=cte) and a variable pressure (p=var) system coefficient kı =0.065 and 0.010 
and a system with volumetric control by a variable capacity pump (Q;=var) 


- determined by simulation 
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Fig. 4. Relation of the overall efficiency n of a constant pressure (p = cte) and a variable pressure (p = var) system with proportional control, with the used 
proportional directional valve coefficient k, = 0.065 and with possible use of a bigger proportional directional valve with k, = 0.010, as well as a system 
with volumetric control by a variable capacity pump (Q, = var), to the load coefficient My, at the value of the cylinder speed coefficient Q = 0.939 
(Vy = 0.380m/s) resulting from the maximum pump capacity Opna: Maximum ng, values of the three considered systems are closer to one another [1]. 


It can be seen in Fig.3 that a 14-fold increase of the cylinder 
speed in the investigated structures causes about 14-fold 
increase of their efficiency. As a comparison, a 14-fold increase 
of the cylinder speed in a Q, = var structure causes about 
2-fold increase of its efficiency (from n = 0.39 at Oy = 0.063 and 
My = 0.875 to n = 0.78 at Oy = 0.875 and My = 0.875). 

Fig. 4 presents efficiency y of a constant pressure (p = cte) 
and a variable pressure (p = var) system with proportional 
control, with the used proportional directional valve coefficient 
k,, = 0.065 and with possible use of a bigger proportional 
directional valve with k,, = 0.010, as well as a system with 
volumetric control by a variable capacity pump (Qp = var), as 
a function of the load coefficient My at the value of the cylinder 
speed coefficient my = 0.939 (vu = 0.380m/s) resulting from 
the maximum pump capacity Qpmax - 

In the maximum cylinder speed range, i.e. with the full use 
of the pump capacity, efficiency values of the p= cte and p = var 
systems with throttling control become close to the efficiency 
of the Q, = var system with volumetric control. 


CONCLUSIONS 


. The maximum achievable values of efficiency of systems 
with proportional (i.e. series throttling) control and of 
a system with volumetric control by a variable capacity 
pump are approximately similar. The compared systems 
were assembled of elements with the same k; coefficients 
of energy losses. 

. By applying a variable pressure (p = var) system, 
a significant increase of the energy efficiency n can be 
achieved with smaller cylinder loads. 

. With small cylinder speed values, the effect of using 
a p = var system is little, mainly due to volumetric losses 
connected with draining the excess liquid to the tank. 

. Optimization of hydrostatic systems means, among other 
aspects, a possibility of foreseeing the behaviour of an 
energy system in various conditions of its operation, as 
a function of speed and load of the hydraulic motor, working 
liquid viscosity, losses in the elements and particularly as 
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an effect of the system structure. The common acceptance 
and use of the objective, experimentally verified methods 
of determining the system energy efficiency, looking at 
the efficiency of entire combined system, can clarify many 
misunderstandings, e.g. those pertaining to the problem 
of maximum efficiency of specific structures. Such as an 
answer to the question: has a system with flow regulator 
different efficiency than a system with common throttling 
valve? 


5. The Author intends to carry out further investigations of 
systems with proportional control, aimed at determining 
the influence of the working liquid (oil) viscosity on the 
energy efficiency. 

NOMENCLATURE 

cte — constant 

fp — throttling slot at the cylinder inlet 

pe — throttling slot at the cylinder outlet 

Fu — hydraulic linear motor (cylinder) load, current force 

required of a linear motor 

Fw — force indicated on the piston of the hydraulic linear motor 

(cylinder) 

Fum — hydraulic linear motor mechanical losses 

Fs» — force of spring in the overflow valve 

k, — coefficient of relative volumetric losses per one shaft 

revolution of fixed capacity pump 

k, — coefficient of relative decrease in pump rotational speed 

k, — coefficient of relative pressure losses (flow resistance) in 

internal pump ducts, at theoretical pump delivery Qp 

k,,; — coefficient of relative mechanical losses in pump, at App, 

=0 

ką) — coefficient of relative increase of mechanical pump 

losses, at increase in pressure in pump working chambers 

k; — coefficient of relative pressure losses (flow resistances) 

in the line joining the pump with throttle control unit, at 
theoretical pump delivery Qp 
k,, — coefficient of relative pressure losses (flow resistances) 
in the line joining the throttle control unit with hydraulic 
motor, at theoretical pump delivery Q,, 

kę, — coefficient of relative pressure losses (flow resistances) in 
hydraulic motor outlet line, at theoretical pump delivery 
Qui 

k,, — coefficient of relative mechanical losses in hydraulic 
motor — cylinder, at a force Fy, = 0 

k,, — coefficient of relative increase of mechanical losses in 

motor — cylinder, at increase of force Fu 

kg — coefficient of relative pressure losses (flow resistances) 

in internal ducts of hydraulic motor, at theoretical pump 
delivery Q,, 

k, — coefficient of relative volumetric losses in hydraulic 

motor 

ky, — coefficient of relative minimum pressure decrease in 

2-way flow control valve, which still ensures the flow 
regulation, or coefficient of relative pressure decrease in 
3-way flow control valve 

kı — coefficient of relative pressure decrease App, in 

directional control valve (servovalve, proportional valve) 
demanded by a maximum throttling section fhgmax for 
receiving flow intensity equal theoretical pump delivery 
Qr 

Mm — hydraulic motor relative load coefficient My= Fy / Fy, 

Po — the reference pressure in the oil reservoir 

Pi — pressure at the cylinder feed proportional valve inlet 

P> — pressure in the outlet conduit from proportional valve to 

the cylinder 

pr  — pressure in the inlet conduit to the proportional valve 

from the cylinder 

P>  — pressure in the outlet conduit from proportional valve to 
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the oil reservoir 
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Pn — nominal (rated) working pressure of hydrostatic 
transmission (hydraulic system) 
Py: — pressure in the inlet conduit to the cylinder 
Pu — pressure in the outlet conduit from the cylinder 
Pmi — pressure in the inlet chamber of the cylinder 
Pmi — pressure in the cylinder discharge chamber 
Pp; — pressure in the pump inlet 
Pp» — pump supplying pressure 
Psp — Operating pressure overflow valve 
Pspo — Opening pressure overflow valve for (Q, = 0) 
Psps operating pressure overflow valve controlled by the 
receiver inlet pressure 
Apco — pressure drop in the inlet conduit to the pump 
Apc, — pressure drop in the inlet conduit to the control unit 
Apc — pressure drop in the line between the control unit and 
cylinder 
Apc; — pressure drop in the outlet conduit from cylinder to the 
proportional valve 
Apc» — pressure drop in the outlet conduit of the cylinder from 
the proportional valve 
Apps — pressure drop in the proportional directional valve 
throttling slot fyy, (at the cylinder inlet) 
Appia — pressure drop in the f,,., proportional valve throttling slot 
(at the cylinder outlet) 
Apu — decrease of pressure (pressure drop) in hydraulic linear 
motor (cylinder) 
Apm; — pressure drop indicated between inlet and outlet chamber 
of the cylinder 
Ap, — increase of pressure in the pump 
AP pp! pressure drop in the inlet channel pump (and the 
distributor, if any) 
AP pp? pressure drop in the pump outlet duct (and the distributor, 
if there is one) 
Q,  — intensity of flow directed through the overflow valve to 
the oil reservoir 
Qu — hydraulic linear motor absorbing capacity, intensity of 
flow to hydraulic linear motor 
Qu — intensity of flow from the hydraulic linear motor 
(cylinder) 
Qp  — pump delivery 
n — energy efficiency 
Su — effective area of the hydraulic linear motor piston in its 
inlet chamber 
Sw  — effective area of the hydraulic linear motor piston in its 
outlet chamber 
SP — overflow valve 
SPS — overflow valve controlled by the receiver inlet pressure 
var — variable 
vu  —_ hydraulic linear motor speed 
wy — hydraulic linear motor speed coefficient — ratio of 
instantaneous speed to the nominal one of a hydraulic 
linear motor — Gy = Vy / Vun 
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ABSTRACT 


In the Faculty of Ocean Engineering and Ship Technology, Gdansk University of Technology, design has 
recently been developed of a small inland ship with hybrid propulsion and supply system. The ship will be 
propelled by a specially designed so called parallel hybrid propulsion system. 

The work was aimed at carrying out the energy efficiency analysis of a hybrid propulsion system operating 
in the electric motor drive mode and at performing the noise pollution measurements. 

The performed investigations have shown that a significant impact on the efficiency and on the acoustic 
emission has the type of belt transmission applied. 


Key words: hybrid ship propulsion; ship propulsion 


INTRODUCTION 


In recent years an unusual development in practically all 
fields of science has been witnessed. In effect, new machines 
and equipment with increasing capabilities are being designed 
and produced and the progress has not left shibuilding out. 

In most of the now constructed ships the conventional 
propulsion with an internal combustion engine is used. Modern 
propulsion units of the type are equipped with complex fittings 
limiting emissions of noxious compounds to the atmosphere. 
Attempts to limit emissions are being made even at the cost 
of increased fuel consumption. Modern internal combustion 
engines are more and more often adapted to use such non- 
typical fuels as vegetable oils or compressed natural gas. 

Specific requirements are imposed on the inland navigation 
ships. Some inland waterways are under particular protection 
as popular recreation regions or due to their exceptional 
landscape values. 

Looking through the shipbuilding technical literature, one 
can often come across the terms like ”no emission”, *”zero 
emission” or ’greenship” in titles of the papers. This testifies 
that the "clean technologies” are now a much more meaningful 
issue than they used to be. 

In the case of a conventional propulsion system, the 
noise pollution” is also a significant problem. Unfortunately, 
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even a perfectly designed engine exhaust system cannot 
fully eliminate noise. Combustion engine is also a source of 
vibration, which may have a negative impact on the condition 
of crew and passengers, particularly during a longer voyage. 
Another problem are the exhaust gases, which with an adverse 
wind can significantly worsen the sailing comfort. 

The electric propulsion has none of the above described 
flaws. It operates silently and generally is not a source of 
vibration. In addition, the most up-to-date solutions have 
a high, 90% exceeding efficiency. However, a problem is the 
motor supply source. The use of a combustion engine-driven 
generating set does not solve the problem, as all the above 
described limitations will then occur. The fuel cell appears 
an advisable solution in the case. Such solution was used on 
a ship already in 1964 [1, 2, 3, 4]. A very interesting ship is 
Alsterwasser launched in 2008. She is equipped with hydrogen- 
supplied fuel cells. 

The same technology has been used also in the 212A 
class submarines built in 1996 by the Thyssen-Krupp Marine 
Systems/HDW corporation to the order of the German Navy. 
They are equipped with the Siemens company hydrogen- 
supplied fuel cells. 

Therefore, the information (given at the beginning of 2013) 
that the Siemens company withdraws from the development 
work on fuel cells was a big surprise. 


Fig. 1. The Alsterwasser ship launched in 2008 in Hamburg, driven 
with a 100 kW electric motor supplied by two 48 kW fuel cells, may take 
a hundred passengers on deck (http://www.zemships.eu) 


Batteries are a commonly used supply source for electric 
motors driving small and medium size inland vessels. 
Unfortunately, the cheapest classical solutions are heavy and 
the modern light batteries (Fig. 2) are very expensive. The 
price of an I kWh module may exceed a thousand Euro and 
it depends on the battery generation, manufacturer, and size 
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Fig. 2. Energy density of the most popular now available battery types 
(based on the Authors’ investigation - data of 2013) 


Therefore, in most of the so far implemented solutions 
the lead batteries have been applied. But this creates serious 
problems. A service ship, e.g. a passenger vessel, is usually 
designed for as intensive operation as possible. It must be 
a source of income to repay, within a shortest possible time, the 
expensive investment of its designing and construction. It is also 
recommended to be relatively fast. Greater speed means greater 
number of trips, besides, the owners’ experience shows that 
passengers get nervous when, even in touristically attractive 
surroundings, the trip takes too long. Greater speed means also 
increased demand for accumulated energy, which translates 
into large volume and great weight of the batteries. Great 
weight of batteries means in turn increased ship draught and 
increased hydrodynamic resistance, i.e. even greater demand 
for energy. Therefore, proper design of the ship hull, with as 
little weight as possible and low hydrodynamic resistance, are 
the basic questions. Very important is also a high-efficiency 
propulsion system. 


As the Authors’ experience indicates, in many cases the 
use of a propulsion system based exclusively on an electric 
motor is not justified. In preparing the ” Design of an ecological 
watercraft for inland waters in Poland” project for realization, 
analysis was carried out of several popular sailing routes, with 
fragments included in the "Nature 2000” programme protection 
system or protected as parts of a national park. It appeared 
that quite long sections of waterways are touristically not 
particularly attractive and there are no navigating speed limits 
on them. Therefore, the project authors came to the conclusion 
that an ideal solution would be a hybrid propulsion system, 
being a combination of a compression ignition engine and an 
electric motor. Such solution allowed considerable reduction 
of the ship weight and the propulsion system flexibility. 


DESIGN OF THE STA.H SHIP 


The designed ship has 19 m length, 5 m width and about 
0.5 m draught and can take on deck from 50 to 70 passengers, 
depending on version. The main dimensions allow it to navigate 
on smaller, shallow Polish rivers, as well as on the route from 
Gdansk to Augustow, where the prototype unit is supposed to 
operate. 

The design of ship is an effect of cooperation between 
engineers from the Faculty of Ocean Engineering and Ship 
Technology, Gdansk University of Technology and artists from 
the Academy of Fine Arts in Gdansk. 

Particular effort was made to reduce the hydrodynamic 
resistance of the ship hull and to minimize draught. Low weight 
was achieved by designing the hull from aluminium. The tests 
of hydrodynamic resistance were performed in parallel for two 
versions — a single hull and a double hull (catamaran) structure. 
Numerical analysis by CFD method and model tests in the 
Faculty of Ocean Engineering and Ship Technology towing 
tank were carried out (Fig. 3). 


Fig. 3. Model of STA.H ship during the towing tank tests; 
a) catamaran b) single-hull version; photo: Mirosław Grygorowicz 
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The tests showed that the catamaran had slightly lower 
hydrodynamic resistance than the single-hull unit. Nevertheless, 
the single-hull version was considered a more promising 
solution due to its smaller draught and smaller overall height 
of the vessel. 

The ship propulsion system has been designed in two 
versions. For fast navigation in the regions with protected 
sections the parallel hybrid diesel-electric solution is proposed. 
In the other version, of a ship designed for navigation in fully 
protected regions, the all-electric propulsion system has been 
chosen. In both solutions the electric propulsion system is 
equipped with a hybrid multi-source supply system consisting 
ofa battery, high power charger, photovoltaic panel module on 
the ship deck and the auxiliary emergency generating set. Both 
solutions have been accepted by the classification society. 


PROPULSION SYSTEM DESIGN 


In the preliminary design stage, the propulsion power and 
the necessary battery weights were assumed from experience, 
as the results of hull hydrodynamic resistance tests were 
not known yet and precise estimation of power demand was 
impossible (Table 1). It has to be remembered, that because 
of the classification society requirements regarding high 
level of safety, two autonomous propulsion systems must be 
provided. In effect, a twin-screw propulsion system must be 
applied, which is otherwise advantageous. The most important 
advantage is smaller diameter of the propulsion screws, 
important for a small-draught inland vessel. Additionally, 
a twin-screw ship has better manoeuvrability. 


Tab. 1. Data assumed in the preliminary design included the weights of 
electric system, propulsion system and supply system 


Brief 


Description foredesign 


Cruising speed with the electric 
motor drive [km/h] 


Power of the electric drive [kW] 


Time of navigating with the electric 
drive [h] 


Battery capacity (with the allowance 
taking into account the loss of 
capacity) [kWh] 


2x 40 


Number of batteries (75kg, 12V, 


230Ah) [pcs.] 3508 


Battery weight [kg] 2 x 1200 


Approximate weight of motors with 


fittings [kg] "A 


After the towing tank tests it appeared that the power needed 
to achieve the required speed is less than the value assumed 
in the preliminary design. However, it has to be remembered 
that in real conditions the impact of wind and of the hull 
overgrowing with water organisms may increase the power 
demand. Also a certain power reserve in the propulsion system 
is needed for manoeuvrability purposes. Finally, a three-phase 
alternating-current electric motor with permanent lantanide 
series magnets was chosen, which, supplied with 48 V voltage 
by means of a dedicated frequency converter, achieves the 
power of 10 kW. It is also important that, if needed, the motor 
can operate as generator driven by the main engine and batteries 
may be charged in this way. 
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The propulsion system designed by the Authors (Fig. 4) 
has some important advantages. It can cooperate with typical 
combustion engines of various power. The designed system 
may be installed aboard a newbuilt ship or an existing solution 
may be modernized. The proposed system may also be used as 
an electric propulsion system without the combustion engine. 
Another advantage is easy assembly and free access to the 
electric motor. In practice, the solution allows easy replacement 
of the drive unit in accordance with customer's demand and 
does not require significant structural changes. 

A complete propulsion system is shown in Fig. 4 below. 
For the system to operate properly, the combustion engine 
must be equipped with a typical reduction gear, with coupling 
allowing to disconnect the engine from the shaft line. Then the 
electromagnetic clutch (4) may be remotely switched on and the 
electric drive activated. It is also worth noting that the solution 
with belt transmission makes it easy to adjust individually 
the transmission ratio depending on the engine type used, by 
exchanging the belt pulleys and belt. 


Fig. 4. Visualization of the STA.H ship hybrid propulsion system: 
1 — combustion engine with reduction gear, 2 — electric motor, 
3 — belt transmission, 4 — controlled clutch 


Similar solutions are now used by established manufacturers. 
Looking through the catalogues, participating in fairs and 
exhibitions, one can hear a lot about advantages of this type 
of drive, particularly about its noiseless operation and high 
energy efficiency. 


INVESTIGATIONS OF A HYBRID SYSTEM 
OPERATING IN THE ELECTRIC MOTOR 
DRIVE MODE 


The propulsion system was assembled in the Faculty of 
Ocean Engineering and Ship Technology laboratory (Fig. 5). 
Combustion engine was not connected to the hybrid propulsion 
system test stand, as the test was aimed at evaluation of the 
propulsion system energy efficiency during the electric motor 
operation. Therefore, concentration was on the electric drive 
together with the belt transmission. The propulsion system was 
loaded with a direct-current generator (pos. 5). 

In order to measure the belt transmission and the electric 
motor efficiency, the controlled clutch (Fig. 4 pos. 4) was 
replaced by a torque meter (Fig. 5 pos. 2). 

The following parameters were recorded during the tests: 
- rotational speed and torque on the motor shaft, 

- torque on the driven shaft, 

- voltage and intensity of the motor input current, 

- noise during the system operation, 

- temperature distribution on the surface of transmission 
subassemblies. 


After completion of the measurements, it was possible 
to calculate the power of losses in the transmission gear and 


Tab. 2. Data of the two investigated drive versions 


First version (I) 


Second version (II) 


Drive type 


Transmission ratio / number 
of wheel teeth 


Motor power / motor 
rotational speed 


Synchronic straight gear belt transmission, 
belt width 30 mm, pitch 8 mm 


Synchronic bevel gear belt transmission, 
belt width 32 mm, pitch 8 mm 


1:3.07 / 80 and 26 


10 kW / 3000 rev/min 


Safety coefficient assumed 
for the drive — service factor 


motor, the propulsion system (motor and frequency converter) 
energy efficiency, the transmission gear efficiency and the 
overall efficiency of the propulsion system. The results are 
presented in chapter 5. 

Details of the two investigated drive versions are given in 
Table 2. 


5 


Fig. 5. Test stand for investigations of the hybrid propulsion system 
operating in the electric drive mode: 1 — electric motor, 2 — torque meter 
installed instead of controlled clutch, 3 — belt transmission, 4 — torque 
meter, 5 — generator, 6 — frequency converter, 7 — resistor 


INVESTIGATION RESULTS 


In designing the drive transmission system, no particular 
difficulties were anticipated and the solution seemed not 
complicated. In the drive first version (I), a high class 
synchronic straight gear belt transmission was used consisting 
of two toothed wheels and one tension roll. In order to be sure 
that the system will operate properly, recommended by the 
manufacturer prefabricated belt pulleys were acquired. The 
small driving wheel had flanged rims protecting against the 
belt slip-off. 

From the beginning of measurements it appeared that 
the drive transmission system operated with excessive noise 
(Fig. 6), the noise intensity reached 87 dB. 

Temperature of the transmission subassemblies was 
measured by means of a thermovision camera (Fig. 7 and 8). 
The measurements were verified by thermocouples (Spl 
and Sp2) installed on the transmission gear. Their position 
was changed in the consecutive measurement series, which 
can be seen in the illustrations. The measurements showed 
that temperature of the full power transmitting drive system 
stabilized after about forty minutes of operation. 


1.5 


During the investigations a new very modern type of 
synchronic bevel gear belt was acquired. Price of such solution 
was almost four times higher than that of the above described 
classic solution. Both solutions were compared. The second (II) 
version had the wheels without protective rims because of the 
specific gear geometry. A comparable efficiency but certainly 
a lower noise level was expected. Results of the measurements 
of temperature distributions (Fig. 8) and efficiency (Fig. 9 to 
11) are presented below. 
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Fig. 6. Diagrams of the measured noise as a function of drive shaft torque 
and rotational speed: a) version I — straight gear belt, 
b) version II — bevel gear belt 
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Fig. 7. Results of temperature field measurements carried out with thermovision camera in the transmission version I; a) photograph of the investigated 
system, b) the system immediately after start-up, c) measurement after 20 minutes, d) measurement in stabilized state after 40 minutes of operation 


21,1 oC 
19,6 


Fig. 8. Results of temperature field measurements carried out with thermovision camera in the transmission version IT; a) photograph of the investigated 
system, b) the system immediately after start-up, c) measurement after 20 minutes, d) measurement in stabilized state after 40 minutes of operation 
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Fig. 9. Diagrams of efficiency as a function of output power for the main shaft speed n = 6.1 rev/s; 
a) transmission with straight gear belt (version I), b) transmission with bevel gear belt (version II) 
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Fig. 10. Diagrams of efficiency as a function of output power for the main shaft speed n = 12.4 rev/s; 
a) transmission with straight gear belt (version I), b) transmission with bevel gear belt (version II) 
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Fig. 11. Diagrams of efficiency as a function of output power for the main shaft speed n = 15.6 rev/s; 
a) transmission with straight gear belt (version I), b) transmission with bevel gear belt (version II) 


DISCUSSION 


The investigation results have shown that the efficiency of 
such belt transmissions increases with the increasing output 
power (Fig. 9 to 11). This is relation typical of that type of 
solutions. The analysis of lost power (Fig. 12 and 13) gives 
a picture of an interesting phenomenon. It appears that the 
power of energy losses does not increase with increasing 
output torque, but slightly decreases. In the Authors’ opinion, 
this may be an effect of stretching the active part of the belt 
and slight slackening of the passive part where the tension roll 
is mounted. In effect, resistance on the tension roll is slightly 
lower. In the case of straight gear belt, the power of losses is 
proportional to the emitted noise (Fig. 12). Such effect cannot 
be seen in the transmission system with bevel gear belt, which 
indicates a different character of the movement resistance. 
In the Authors’ opinion, the main source of losses and of the 
acoustic emissions in the case of straight gear belt is throttling 
the tooth space by consecutive belt teeth and making the 
compressed air come out sidewise. This effect can be limited by 
removing one flanged rim of the belt pulley or making holes in 
the rims to let the compressed air out. In the investigated case, 


removing one rim was not a promising method as the belt was 
several millimetres narrower than the pulley. 

Worth considering is also a conception of not 
overdimensioning the drive system by assuming a high service 
factor. In designing the analysed system (after consultation with 
the manufacturer), the value 1.5 of service factor was assumed. 
In this way it will be possible to extend the transmission gear 
operating field in the range where the efficiency exceeds 90% 
(Fig. 14). However, this can be achieved at the expense of the 
drive durability. In the case of a typical inland vessel, which 
usually has a long winter break in navigation, it does not seem 
to be a problem as the belt can then be renewed. 

In the analysed case (Fig. 14), decreasing the assumed 
service factor from 1.5 to 1 would cause extending the operating 
field where efficiency exceeds 90%. 

The investigations have shown that the bevel gear belt has 
higher efficiency than the straight gear belt by about 5% in 
the whole range of rotational speed and torque where the tests 
were carried out. A considerable surprise was noisy operation 
of the straight gear belt drive system (Fig. 6 and 12), which 
is surely an effect of the high rotational speed of small wheel 
reaching 3000 rev/min. 
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Fig. 12. Diagrams of lost power and measured noise of the drive system with straight gear belt (version I) as a function of the drive shaft torque 
and rotational speed: a) main shaft rotational speed 15.6 rev/s, b) main shaft rotational speed 12.4 rev/s 
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Fig. 13. Diagrams of lost power and measured noise of the drive system with bevel gear belt (version II) as a function of the drive shaft torque 
and rotational speed: a) main shaft rotational speed 15.6 rev/s, b) main shaft rotational speed 12.4 rev/s 
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Fig. 14. Diagram of efficiency as a function of output power 
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SUMMARY AND CONCLUSIONS 


The energy efficiency ofa battery-supplied ship propulsion 
system is of primary importance. Batteries have limited capacity 
and the accumulated energy must be utilized economically. 
Cases have been known where crews accustomed for years to 
operate ships with conventional propulsion could not cope with 
a battery or photovoltaic cell-supplied electric drive of ships. 

In analysing the propulsion system overall efficiency it 
has to be remembered, that it is a product of the efficiencies 
of motor, transmission gear, bearings and screw propeller. If 
in a ship a high efficiency screw propeller with maximum 
efficiency of 65% were used, then in optimum conditions the 
overall propulsion system efficiency would not exceed 55 %. 

The tested 11 kW electric motor reached efficiency 
exceeding 90% in the output power range above 5.5 kW. Power 
of only about 10 kW was obtained during the tests. The motor 
did not achieve the catalogue rotational speed and torque. The 
Authors had experienced similar situations in the past, when 
the drive unit parameters differed from the catalogue values by 
more than 10 %. In some cases the operating parameters could 
be improved by reprogramming the frequency converter. 

The investigations will be continued in order to compare 
properties of other two types of belt transmissions. Transmission 
with five wedge belts and transmission with a multi-wedge belt 
have already been designed. 
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ABSTRACT 


With the substantial rising of international oil price and global warming on the rise, how to reduce 
operational fuel consumption and decrease air pollution has become one of the pursued goals of green 
ship. Ship route planning is an indispensible part of the ship navigation process, especially in transoceanic 
crossing ship routing. The soundness of ship routing not only affects the safety of ship navigation but also 
the operation economy and environmental protection. This research is based on the platform of Electronic 
Chart Display and Information System (ECDIS), and founded on Ant Colony Algorithm (ACA) combined 
with the concept of Genetic Algorithm (GA), to model living organisms optimization behaviour to perform 
efficient ship route planning in transoceanic crossing. Besides the realization of route planning automation, 
ship routing will achieve the goal of optimum carbon dioxide reduction and energy conservation, and 
provide reference for route planning decision. 


Key words: ship routing; Ant Colony Algorithm; GIS; weather routing 


1. INTRODUCTION 


Almost 90% of the world’s trade is carried by ship and, for 
the vast majority of this trade there is little or no alternative 
to transport by ship. Statistical data indicate that the freight 
carried by global commercial shipping has reached 8.919 
million tons, and that this large volume has caused substantial 
negative environmental impacts, even though ships are the 
most environmentally friendly form of transportation. Still, it 
is estimated that 4-5 per cent of global CO, emissions come 
from international shipping [6]. According to statistics, there are 
about 52,000 registered civil ships of more than 100t capacity, 
annually consuming 500 million tons of fuel, producing 1.6 
billion tons of carbon dioxide and 30 million tons of sulfur 
dioxide, causing greater global warming and air pollution. 
Therefore, an additional option for reducing air pollution is to 
decrease fuel consumption during transit [3]. From a climatic 
and environmental perspective sustainability has become an 
issue of the highest priority, and this is an imperative that 
cannot and should not be ignored. Within the foreseeable 
future shipping will still be dependent on fossil fuels. We 
must develop energy efficient and environmentally friendly 
technologies encompassing shipboard energy production, ship 
propulsion and ship operation. In recent years the concepts of 
the “green ship” and “green shipping” have been raised in the 
international shipping industry. The green ship concept covers 
the full ship life cycle, and provides total solution to pollution 
prevention, noxious gas emission, and energy conservation 
problems associated with the ship’s design, production, 
operation, and scrapping. Therefore, shipping management 
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should take care not only of economic efficiency, but also 
environment protection to prevent pollution. The coordination 
of efficiency and environment allows for environmentally 
sustainable development. In order to achieve this goal, ships 
already in operation can be retrofitted with energy-conserving 
equipment but the most direct and economical method is to 
plan an optimum route. Choosing the optimum route, adjusting 
the speed and avoiding adverse weather can not only reduce 
fuel consumption significantly but also improve navigational 
safety. 
The efficiency of transoceanic ship route planning is 
generally determined by the following three aspects [26]: 
- the accuracy of the prediction of the ship’s hydrodynamic 
behaviour under different weather conditions; 
- the accuracy of the weather forecast; 
- the capability and practicability of the optimization 
algorithm. 


This paper deals with navigation technology and focuses on 
the third above mentioned aspect, i.e. ship route optimization. 
Because optimum route determination is a multi-criterion 
non-linear planning issue and encompasses many constraints, 
a balance between navigational safety and economic efficiency 
should be found [18]. In other words, the applied planning 
framework not only has to take care of navigational safety- 
related risk evaluation and deviation but also avoid total cost 
increase due to deviation measures. In order to achieve this 
optimization mission we made use of the new Ant Colony 
Algorithm (ACA) developed in the evolutionary computation 
field to model and simulate how ants choose foraging paths. 


We modified the ACA for marine navigation application, and 

integrated the Genetic Algorithm (GA) concept to improve 

its practical effectiveness. The geometric computation and 

information search function is provided by the Geographic 

Information System (GIS) - based platform ECDIS along with 

weather information for the background data required for ACA 

calculation. Through this research design we hoped to achieve 

optimum ship routing and use a more objective, quicker and 

more accurate model in order to: 

- automatically generate a more rational recommendation route 
that gives consideration to both safety and economy, 

- lighten the navigator’s workload, 

- provide a reference for route planning decision making, 

- collaborate in the green shipping and environment 
sustainability effort. 


2. PREVIOUS WORK 


Since research in transoceanic crossing weather routing 
typically studies broader geographic domains, usually across 
the Atlantic or Pacific, it mainly considers mid- and short-term 
weather effects, especially the effects of wave direction and 
height on a ship’s cross-ocean navigation speed, to achieve the 
goals of minimum sailing time and fuel consumption. After the 
appearance of an innovative paper by Hanssen and James [8], 
several methods have been proposed to solve the problem of 
minimizing sailing time. The major methods include: 

1. Isochrone method. An Isochrone Line (Time-Front) is 
defined as the outer boundary of the attainable region from 
the departure point that is reached after a certain time. The 
isochrone method proposed by Hanssen and James [8] was 
used for along time because it offered an easy manual 
calculation method. Hagiwara [7] devised a modified 
isochrone method suitable for computer programming. 
However this method can generate so-called “isochrones 
loops”, and cannot ensure that the route will not cross 
land. Hence to counter these problems, Szlapczynska and 
Smierzchalski [19] made improvements to the isochrone 
method, consisting in generating the initial set of candidate 
routes necessary for evolutionary computation to enhance 
the effectiveness of initial random population generation 
in evolutionary computation. However this method utilizes 
bitmap-based algorithms and is not applicable to vector- 
based ECDIS. 

2. Calculus of variation. In Bijlsma publications [1, 2, 3] 
application of the calculus of variation was extended to 
the minimization of fuel consumption. Kosmas et al. [11] 
used calculus of variation and a wave model to compute the 
optimum route and applied simulated annealing techniques, 
a type of optimization method, to accelerate algorithm 
convergence. This approach is limited to theoretical 
analysis, time-consuming in terms of calculation, and less 
applicable in practice [24]. Khalilov [10] states that the 
variation calculus method has no optimality proof. 

3. Dynamic programming. Dynamic programming uses 
a grid system which divides a possible sailing region into 
several cells. Each crossing point of a cell boundary is 
a waypoint candidate. The dynamic programming solution 
algorithm seeks a ship’s trajectory that is composed of 
the positions (X, Y T),, k= 1, 2,...,N, and the Control C, 
with the initial time 7, and the initial ship position X [12]. 
Motte et al. [14, 15, 16, 17] conducted a series of detailed 
investigations on this issue, including grid mesh size 
setting, system development and design. In recent years, 
Wei et al. [25], Zhou and Chen [27] and Tang et al. [20] 
have also used dynamic programming to conduct automatic 


ship routing to reduce the workload of navigators. In 
Bijlsma paper [2] the relation between optimum control 
theory and dynamic programming was elucidated for 
the case of minimal fuel routing. The development of 
this method has become more mature but its weaknesses 
include the large quantity of nodes, large amount of space 
required for storage, inflexibility, and slow computational 
speed. 

4. Network Model. A network is a special form of graph that 
consists of nodes and arcs. The network has numerical 
parameters such as length, cost and transit time. Lee et al. 
[12], through network modelling, used a depth-first search 
algorithm to combine weather modelling to establish 
a web-based optimum ship routing system. Montes [13] 
modelled the Optimum Track Ship Routing (OTSR) for 
U.S. Navy warships using a network graph of the Western 
Pacific Ocean. A binary heap version of Dijkstra’s algorithm 
determines the optimum route for given model-generated 
wind and sea inputs. This method is similar to the dynamic 
programming method, but computation is inefficient 
because of a large quantity of nodes. 


Problems such as difficult computerization, complicated 
and inefficient calculations, and easily getting trapped into local 
optima, are found in the above-mentioned traditional solution 
methods. With the development of artificial intelligence and 
the maturing of evolutionary computation technology, more 
and more researchers investigate model organism optimization 
by using evolutionary computation technology (for example: 
genetic algorithms (GA), ant colony algorithms (ACA) and 
particle swarm optimization) and apply it to route planning 
solutions. Tsou [22, 23], Smierzchalski and Michalewicz, 
[18] and Ito et al. [9] used the GA in the planning of collision 
avoidance and coastal navigation paths and gained very good 
results. The present research used the ACA as the research 
method, emulating ants’ forage path finding behaviour, and 
applied it to search an optimum route. The two problems are 
similar with regard to their nature; ships are considered as 
ants, optimum navigation routes as forage paths, and weather 
conditions as obstacles. We hope that through this conversion 
approach the ACA calculation efficiency and global optimum 
search capability can be enhanced in achieving optimum ship 
routing. 


3. OPTIMUM SHIP ROUTING AND 
WEATHER ROUTING 


The ship weather routing develops an optimum track 
for ocean voyages based on ship’s individual characteristics 
for a particular transit, and forecasts of weather and sea 
conditions. The advent of extended range forecasting and the 
development of selective climatology, along with powerful 
computer modelling techniques, have made ship routing 
systems possible [4]. 


3.1. The principle of optimum ship routing 


Today, in the evaluation of optimum transoceanic routes, 
the primary concern, besides safety, is usually given to 
navigation time (or navigation speed) rather than to distance. 
In synthesized considerations of safety and energy conservation 
the shortest distance between two ports is not necessarily an 
optimum route. If there are many routes to choose from, under 
the condition that the main engine power outputs are equal, the 
route that takes the least amount of time is the optimum route. 
Shorter navigation time in itself implies less fuel consumption. 
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Therefore, the optimum route must consider at least the 

following conditions: 

- safety; 

- avoiding severe weather and sea conditions as much as 
possible; 

- natural factors such as weather, sea conditions and ocean 
currents, to increase navigational speed as much as 
possible. 


In order to find the optimum route, when choosing and 
planning the navigation route and making corrections, the 
following principles must always be met: 

- firstly, to ensure safety, avoid severe weather and sea 
conditions as much as possible; 

- secondly, if, under the circumstances, there is no other 
option but to pass through a severe sea environment, pass 
in the least time possible; 

- thirdly, the design of navigation routes always aims to 
shorten navigation time and arrive ahead of the ETA 
(Estimated Time of Arrival); 

- fourthly, as much as possible make use of natural conditions 
which help to increase navigation speed, such as down 
winds and down currents, and as much as possible avoid 
natural factors which work against speed of navigation, 
such as head seas. 


The impact of wind and waves on ship speed is non-linear. 
In circumstances where other conditions are unchanged, their 
effect on ship’s speed, main engine power outputs, and fuel 
consumption cannot be represented by a simple algebraic 
expression. Therefore, in terms of ship’s energy conservation, 
if, upon planned corrections, navigation overall makes use of 
natural conditions that work to its advantage, and avoids natural 
conditions that work against it, total navigation time usually 
can be decreased, although the total navigation distance may 
be lengthened, because increased navigation speed can usually 
compensate for an increase in navigation distance. 


3.2. Computation of optimum route 


Based on the above-mentioned principles of optimum 
route planning, the cost of the optimum route is a function of 
the extent of obstruction and the distance between start and 
end points. Therefore, solving an optimum ship route problem 
may be seen as finding the route with the least cost among the 
set of candidate routes, while taking weather conditions into 
consideration. The mentioned cost is not limited to distance 
only but may include navigation time, fuel consumption or 
a combination of factors. The decision variables that should be 
controlled are the manoeuvring of course and speed. The total 
computational cost (I) is calculated as the sum of the costs over 
each route leg (S). It is based on the ship’s position (P), changes 
in the manoeuvring control (C) (i.e., speed and course) and the 
time duration (t). With these cost components, the ship routing 
problem can be described by the following formula: 


to minimize: 
I =|f (P(s), C(s), t(s))ds (1) 


where: 
f(P, C, t) — function of voyage cost for position, manoeuvring 
control and time; and Pe R, Ce Cy 


Let R be a possible sailing region. The manoeuvring 
control variable C, is the limit of the ship’s output speed and 
the limit of course changing range. The cost per unit distance 
ds is measured by fuel consumption and degree of ship safety 
under sailing. 
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4. ANT COLONY ALGORITHMS AND SHIP 
ROUTE COMPUTATION 


4.1. The Ant Colony Algorithm principle 


The ant colony algorithm is a form of evolutionary 
computation that models the behaviour of real ants in search of 
food. The algorithm was proposed by Dorigo [5] and has been 
successfully used to solve many real-life problems, such as the 
travelling salesman problem (TSP). In the natural world, in 
ants’ search of food, among many different possible routes from 
the nest to the food source, different ants will at first choose 
different paths but in the end nearly all ants will discover the 
same shortest path. This is because ants’ search for the shortest 
path is an interactive process. They leave a certain amount of 
secretion (pheromones) on the path they pass, can sense the 
presence and intensity of the secretions and move toward the 
direction where the substance has a higher concentration. The 
amount of secretion increases with the increase in number of 
ants that pass by, and decreases in a certain function of time 
that passes by. The reason is that more ants pass by the shortest 
path and the speed of pheromone accumulation is higher 
than that of other paths. Through the continuous exchange 
of pheromones the optimum route between the nest and the 
food source is found in the end. The ACA is an optimization 
decision method designed according to this ability of ants in 
the natural world. One sets a certain number of artificial ants 
(ships) and follows a set of search rules to find the optimum 
route comprehensively. The objective of weather routing is to 
find the safest and most economical navigation route based on 
the weather and manoeuvring control conditions between the 
departure and destination points. The conditions and objectives 
are basically the same as those of the ACA. 


4.2. Division of the marine grid system 


In order to plan the necessary deviational and directional 
change in the navigation route, maritime space must first be 
divided into grids. The grids form the basis for waypoint 
candidates, and the navigation route structure consists of 
the waypoints (Fig. 1). Motte [16] has conducted detailed 
comparison and analysis regarding grid mesh size setting, 
and suggested that the grid mesh size should be four degrees 
longitude by a half degree latitude. 

However the experimental measures are based on the North 
Atlantic Ocean from Western Europe to North America and may 
not be applicable to other marine areas. Additionally, given the 
case of a near north-south travel direction, the great circle route, 
using four degrees longitude as the grid’s horizontal width, is 
not always appropriate. Therefore, in this study distance is used 
as the basis for the setting of grid size, hence flexibility for users 
themselves to change settings is saved. Essentially, navigation 
speed should be referred to, so that grid size corresponding to 
about the distance of one day’s navigation, can be kept. The 
procedure to make a grid is as follows: 

- first establish a great circle route between the starting point 
(L,, 4,) and end point (L,, à), to be the basis for route 
planning. The Great Circle route distance (D) and the initial 
course (C) can be obtained from the following formula: 


CosD = +SinL,-SinL, + CosL ,-CosL,-CosDLo (2) 
+ CosL, : TanL , — SinL , : CosDLo 


Coit = - (3) 
SinDLo 

where: 

DLo — the difference of longitude. 


When latitude L,, L,have same sign, then the signs of SinL, 
SinL, and CosL, TanL, are (+); otherwise they are (-). 

- then set a distance X, and along the great circle route set 
segmentation points (Lx, Ax) at intervals of X. The location 
(Lx, Ax) ofa segmentation point from a departure point with 
a certain distance D, and initial course C is calculated by 
the formula: 


L,=Sin '(SinL,- CosD, + CosL,-SinD,-CosC) (4) 
X 1 x 1 X 
CosD , —SinL, -SinL , 
A IO] 
Cosh, “Cosh, 


- then, using each segmentation point as a basis, extend 
the longitude line to the north and south and pick a set 
distance Y to establish each segmentation point, so that 
the coordinates have the same longitude and different 
latitude; 

- then connect these segmentation points in the north-south, 
east-west direction to form a grid system. 


Ay = à + Cos" ( 


Finally set up the coordinate and node matrix, and the 
route planning will be established on this grid system. With 
advancements in communication and positioning technology 
the grid system and navigational route can be dynamically 
adjusted continuously, according to the latest weather 
conditions and the latest positioning point and end point on 
the great circle route. 


4.3. Implementation of Ant Colony Algorithm 


First, establish a grid system and its node matrix in the 
navigation area, and provide each node with the appropriate 
initial value in the matrix to form the initial node information 
concentration value matrix. Then, put all the ants on the starting 
point to allow them to move simultaneously toward direction 
of the aim, and to reach the destination point in the end (Fig. 
1). In the process every ant uses the state transition rule and 
undergoes a selection process when gathering at the next node. 
Assuming that the time it takes for an ant to travel from node 
a in the node set i to any node b in the next node set i+1, is the 
same, and distance is not a factor here, all the ants will reach 
the destination simultaneously and complete a circle. After all 
the ants have reached the destination, calculate the objective 


function value (Eq. 1) based on the routes discovered by 
each ant. Each node is also updated with overall information 
quantity, and nodes that are not passed by are treated as 
information quantity evaporation. Repeat this process until the 
optimum route is discovered. The detailed explanation of the 
process is as follows: 


4.3.1. Adjustment of information quantity 


The number of ants is assumed m. And, every ant chooses 
the next route according to the corresponding probability based 
on the information quantity on the route, instead of choosing 
the route it has previously completed. Upon completing a route 
circle, it produces a corresponding quantity of information 
concentration, based on the total distance of the route, and 
updates the information concentration on the route just 
passed. 

With T, (t) representing the information concentration on 
route (a, b) at t-moment, the information concentration on this 
route at t+] moment is: 


m 
k 
Tap (t +1) = P -Tap (t) + 2 ATi (6) 
k=1 
where: 
p  — a constant which indicates that the residual factor of 


information quantity is between 0 and 1, and (1 — p) 
represents the information evaporation coefficient 

Atk, — the increased information quantity of ant k on route 
(a, b) between t and t+1 moment 


_Q If the k-th ant passes 
, |d, through a, bin this cycle 
M >| (7) 


0 else 


Eq. (7) indicates that the increased information concentration 
is related to the cost of the passed route, and: 

Q is a constant which is the total information amount 
produced by each ant; d,, is the distance between node a and 
b, which can be calculated through Mercator Sailing, and 
then expressed in terms of other costs, such as time, fuel 
consumption etc. 


: Deśtination 


Fig. 1. Grid system and Great Circle Route 
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4.3.2. Waypoint selection principle 


The probability of the event that the ant k which is at node 
a will select node b as its next destination at t moment, is: 


Tap (©) = Nab 
Litt) Mh if bh EN 


PE (W = 4 DEN a) k(a) (8) 
0 else 

where: 

Nka) the set of nodes that are yet to be visited by the k* ant 
at node a 

Na» — heuristic information transferred from node a to node 
b. The information is obtained through the problem 
which has to be solved. In the problem the reciprocal 
of the cost of the points a and b, is n,,= 1/d,, 

a  — the relevant importance degree of residual information 
element concentration 

B — the relevant importance degree of the expected 


value 
4.3.3. Calculation procedure 


The steps for using the ACA to determine the optimum 
route are listed below: 

1. Initialize the information quantity on all the nodes on the 
navigation area, to form the initial information quantity 
matrix. 

2. m ants in number are ready to start at the starting point A. 

3. Every ant selects the next node on the matrix based on the 
state transition rule of Eq. 8, and reaches the destination in 
the end, forming a feasible route. 

4. Calculate the objective function of every ant’s feasible route 
based on Eq. 1, and keep the optimum route. 

5. According to the objective function and the information 
quantity adjustment principle of Eq. 6, adjust the information 
quantity of every point. 

6. Check if iteration conditions are met (if the established 
iteration number or minimum objective function value 
is reached). If the conditions are met then the search is 
complete. If not, repeat the procedure starting from the step 
(2) until the iteration conditions are met. 


4.4. Calculation of navigation cost 


The calculation of navigation cost must take environmental 
conditions and ship performance conditions into account. To 
simplify the numerical experiment the amount of speed loss is 
calculated only by using the effect of waves. 


4.4.1. Ship response and performance 


Ships travelling in the ocean are affected by wind and waves. 
As a result the actual speed of the ship in waves will certainly 
be lower than that in calm water. The degree of speed loss will 
affect the calculation of the speed-weather conditions as well 
as computed navigation time or fuel consumption costs for the 
ship route [14]. Speed loss is the most important among the 
various factors involved in ship route calculations. It is present 
in many route calculation formulae and has significant influence 
on the precision of vessel positioning and the feasibility of 
route optimization by weather routing. In general, there are two 
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approaches to studying ship speed loss. The first approach is to 
conduct seakeeping tank and wind tunnel model tests to simulate 
air and water resistance on the ship body and to derive theoretical 
expressions to quantify the effect of wave resistance on the ship. 
The second approach involves the use of mathematical statistics 
or neural network regression to obtain empirical formula for 
calculating ship speed loss due to waves, by means of a large 
database of weather observation information (e.g. weather 
records or logbook) and information regarding ship conditions 
at the time. The formulae for speed loss obtained from the two 
methods can be used to produce ship performance curves or 
speed loss curves, as shown in Fig. 2. The effect of wave direction 
on speed loss is different for head seas, following seas, and 
beam seas. The relationship between ship orientation and wave 
directions is demonstrated in Fig. 3. Based on ship performance 
curves in combination with sea condition forecast data such as 
wave height, wave direction, or wind speed, the impact of ship 
speed loss on the total navigational costs can be calculated. 
Empirical formulae for ship performance curves have been the 
subject of a number of studies in the past, where different factors 
were considered. The main factors affecting ship speed loss are 
wind and waves, and wave size is closely related to wind speed. 
For simplicity of formulation and because the present study is 
focused on route optimization algorithms, the empirical formula 
recommended by Soviet Central Maritime Research Institute 
to calculate ship speed loss has been selected as the basis for 
comparison. The included main factors are wave height, wave 


0 20 40 60 80 


Significant wave height (feet) 
Fig. 2. James ship performance curves (Motte et al., 1988) 
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Beam seas 


Own Ship 


Following seas 


Fig. 3. Ship orientation and wave direction 


direction, and ship performance coefficient. If necessary, this can 
be substituted in the future by more sophisticated and accurate 
formulae. However, the assumed formula is sufficient for the 
purpose of validation of this research results. The formula is 
expressed as follows: 


V=V-0(0.745h-0.257qh)(1.0—1.35x10 DV,) ©) 


where: 

— the actual speed in the sea 

o — the speed in calm water 

the wave height 

— the angle between ship heading and wave direction 
— the actual displacement of the ship (tons) 


Ue F<4< 
| 


This formula is applicable to ships with displacement 
ranging from 5000 to 25000 tons, and speeds between 9 and 
20 knots. Research has shown that when the wave height 
is between 0 and 5 meters, the difference between speed 
calculated by using this formula and the actual speed is less 
than I nautical mile. However, the calculated speed becomes 
slightly higher than the actual speed when the wave height is 
between 5 and 8 meters. 


4.4.2. Environmental data processing 


Environmental data include geographic, navigational, 
and weather data, i.e. the basic background information on 
navigation areas, and the primary basis for cost calculation. 
By means of intersection calculation and measurement of 
navigation distance and environmental information, different 
cost data can be obtained. The processing and calculation 
method is as follows: Geographic data. Land and ocean 
distribution and danger zone (e.g. shallow waters) distributions 
are stored as GIS’s polygon vector data. With regard to land 
area, overlap analysis is conducted on the marine grid system's 
node and land distribution. Nodes in land areas are initialized 
with the value of 0, while other nodes are set at a constant 
value. Thus, ants will not enter the nodes in the land area and 
will only search in the marine area. 

- Weather data. With advances in information and 
communication technology in recent years, onboard internet 
connection by the INMARSAT ship station is often available 
and ample weather information is then easily accessible. 
Because ocean-going vessels can be only connected to 
the internet by satellites, which is relatively costly, only 
few ships use this method to gather weather information. 
However for the optimum route planning function in the 
Electronic Chart Display and Information System (ECDIS) 
the application of this method is an inevitable trend. In 
this study, the GRIB data downloaded from the internet is 
converted to a GIS format in the system and used for the 
route optimization. 


The GRIB hydro-meteorological data is the standard defined 
by the World Meteorological Organization, and many service 
stations in the United States and Europe publish meteorological 
data in this format. In general, forecasts for the next 5 days 
are available, and they are updated every 6, 12 or 24 hours. 
The spatial resolution can reach 0.5” x 0.57. This format offers 
great flexibility in that wind speed, wind direction, air pressure, 
seawater temperature, wave height, wave direction, and some 
other meteorological information can be selectively included. 
Such information adequately meets the requirements of weather 
routing. The information can be downloaded directly from the 
internet to the computer onboard for further processing, saving 
time and improving accuracy. Because the loss of ship speed is 


mainly caused by waves and Eq. (9) is taken into consideration, 
in this study wave direction and wave height are extracted from 
the GRIB data for subsequent calculations. The wind effect was 
ignored because its relationship with waves is linear. 


4.4.3. Calculation method 


Wave direction and wave height grid data extracted from 
the GRIB are subsequently converted to cell layer data in the 
GIS to facilitate the integration with GIS computation. The 
weather cells are different from the nodes which constitute 
navigational waypoints in the navigation grid system. Next, the 
use of GIS's spatial analysis function is made to calculate the 
distance across all the cells that every candidate route passes, 
and find the wave height and wave direction data of the passed 
cells. By using ship speed loss formula the speed loss quantity 
can be computed. Finally, based on the distance across the cells, 
time/fuel consumption cost can be calculated. If navigation 
time is considered the primary cost, the computation formula 
is as follows: 


d; 
t= (10) 
| y-w, 
T,= >); +4, xF) (11) 
jeJ 
where 
t — the required time for candidate route to pass the j" cell 
d; — the distance for candidate route to pass the j* cell 
v — normal speed without weather disturbance 
w, — the quantity of speed loss due to waves in the j* cell 
T, — total sailing cost 
F; — additional unit distance cost factor in the j* cell 


4.5. Improvements in computational performance 


4.5.1. Limitations of maximum course deviation and 
maximum route segment length 


Because a great number of nodes in the vast ocean grid 
system is not feasible in relation to navigation planning without 
proper heuristic calculation, the process is very time-consuming 
and much calculation capacity will be wasted. Because of 
limitations of maximum course deviation and maximum route 
segment length much unnecessary route calculation results may 
be spared. The course deviation angle is limited to + 60° from 
the original course and the maximum route segment length 
cannot exceed that which can be covered in two days with the 
current speed. The selecting of appropriate values is a topic of 
further discussion. 


4.5.2. Determination of critical speed 


While performing ship routing calculations we cannot seek 
only the shortest navigation time and must also consider a ship’s 
critical speed (allowed maximum speed) while navigating in 
storm conditions. Thus the calculation result can produce the 
optimum route on the premise of safety. The calculation of critical 
speed is mainly related to wave height and wave direction. In this 
study the following formula [24] is used to impose maximum 
speed limits on a ship’s critical speed in waves. 


Vin = e Bia- | r(q) (12) 


imit 
where: 

u(q) = 12.0 + 1.4x10-4q23 
r(q) = 7.0 + 4.0x10-4q23 
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stands for the relative wave direction 
stands for wave height 


q = 
h = 


The formula reflects ship’s critical speed under different 
wave conditions to ensure navigational safety. 


4.5.3. Limit on information quantity 


When the scale of the problem increases, the information 
quantity of the grid point which has never been searched, 
will decrease near to 0, because of the information quantity 
evaporation coefficient (1-p), decreasing the algorithm’s 
global search capability. If (1-p) becomes too small then 
the probability of a solution which has been searched, being 
selected, may become too large and will also affect the 
algorithm’s global search capability [21]. Although increasing 
(1-p) may enhance an algorithm’s global search capability 
but also decrease the algorithm’s rate of convergence. Thus, 
a range of information quantity (Tins Tmax), is established, and 
if the information quantity exceeds this range the following 
corresponding limits are imposed: 


Tmin Tap (t + 1) < Tmin 
T4(t+1) =5 Tmax Tbt] St. (3) 
Ta (t +1) other 


4.5.4. The addition of crossover operations to the genetic 
algorithm 


The crossover operation concept is added to the GA in 
this study. The crossover operation is processed through the 
searched optimum route and another ant's random route every 
time. If the two ants pass by the same grid point (waypoint), 
then a crossover operation is conducted between the two 
routes (excluding the starting and end points) based on the 
same grid point, namely, by combining the first half of route A 
with the second half of route B to form a new route, and also 
the second half of route A and the first half of route B to form 
another new route. If the two ants do not pass the same grid 
point then crossover operations are not conducted. If, after the 
crossover operation, a route better than the current optimum 
route is obtained then the optimum route information is updated, 
otherwise, another crossover operation is performed. The 
adding of crossover operations increases the ant's ability to find 
a better route in a search, as well as variety of solutions. 


4.5.5. The addition of mutation operators to genetic 
algorithms 


In order to avoid local optima, in this study the concept 
of mutation operation was added to the GA through setting 
a certain probability value randomly selecting a grid point 
from the current optimum route and exchanging it with a not- 
yet-passed random grid point from the optimum route, to form 
a new route. If the new route has a lower cost than that of the 
original optimum route then the original optimum route is 
replaced with this route. 


5. EXPERIMENT RESEARCH RESULTS 
5.1. Condition settings of the experiment 
This experiment research simulated the navigation route 


from Yokohama to San Francisco, with (34°40' N, 140°E) as 
the starting point coordinates, and (37°45' N, 122°W) as the 
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end point coordinates. The following conditions were applied 

to the simulation. The time spent in the algorithm is verified 

in the results. 

- Weather conditions: the presumed weather conditions in the 
research had been previously ascertained. The weather cell 
size of wave height - wave direction weather data was set 
at latitude 0.5°I longitude 0.5°. The GRIB data used was 
for the northern Pacific Ocean on 2011/11/25. 

- Conditions of course deviation and route segment length 
restriction: it was assumed that the ship was located at the 
X node (waypoint) in the node set I, and the ship’s next 
feasible node (waypoint) was confined within the range of 
X + 3 node (waypoint) in the node set I + 1. This approach 
was used to provide the conditions of course deviation and 
route segment length restriction. 

- Ship performance was simulated under the output power 
required to maintain the navigation speed of 18 knots in 
calm water. The speed loss curves, as in Eq. 9, were matched 
with wave height and wave direction data to simulate ship 
speed. In the simulation the ship displacement of 18000 
tons was assumed. 

- Grid system size: as shown in Fig. 1, the east-west direction 
interval distance along the great circle route was set at 
300 nautical miles, with a total of 14 points. The distance 
along the longitude north-south direction line was set at 
60 nautical miles (1°), with a total of 44 points. Thus, 
a grid system consisted of 616 waypoints was established. 
The route was generated through the ACA search; it was 
generated by testing and combining the waypoints. 

- ACA parameter setting: a = 3, B = 2, p= 0.5, Q = 1000, m 
(number of ants) = 100, N (number of cycles) = 100. 

- Experiment environment: to perform the simulations PC 
of 4 GB memory and an Intel Core 2 Duo CPU 2.2 GHz 
processor, was used. Visual Basic 2010 was used as the 
programming language and the geographic information 
feature provided by ArcGIS Engine’s COM object was used 
to perform spatial analysis and geometry computation. 


5.2. Discussion of results 


The shipping routes were assessed on the basis of the 
total navigation time. From the experiments it was found that 
appropriate choice of the number of ants and the parameters 
can improve the efficiency of the algorithm. When the number 
of ants is small (ant no. 10), the ACA does not show clear signs 
of convergence, and an optimum solution is hard to obtain (as 
shown in Fig. 4). This is also the case when the ant number is 
larger (ant no. = 20) and parameters a and p are of low value 
(both < 3). When a and B are larger (both > 4), convergence 
is more evident, but a large number of iterations is required. 
When ant no. > 50 and a, B both > 4, a better convergence is 
achieved; however, the solution fell into a local optimum on 
many occasions (see Fig. 5). 

It was found that when the number of ants is relatively 
large and the parameters are such that a < B (a= 2, B= 3), 
better convergence was achieved through the ACA and a near- 
optimum solution was more readily obtained. Therefore the 
parameters for deriving the optimum shipping route are set 
as follows: 


a=2,B=3,p=0.5, m= 200 and N = 100 


In addition, in order to validate the results the simulated 
results from east bound and west bound navigation were 
compared with those from navigation following the Great 
Circle Route. The performance of the simulation program 
under different weather conditions was observed. The simulated 


Ant no.=10 


Hours 


Iterations 
Fig. 4. Poor convergence for a small number of ants. It is difficult to obtain 
a near-optimum solution. Notation: Ant no. stands for number of ants 


(Fig. 5) 
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Fig. 5. When a = 3 and B = 4, the solution achieved through ACA is often 
the local optimum. Notation: Ant no. stands for number of ants 


results are shown in Fig. 6 where the purple route is the Great 
Circle Route, the black route is the optimum route produced 
by the ACA, and the yellow routes are the alternative routes 
generated in the search process. During the east bound 
navigation, due to the influence of strong winds and stormy 
seas in the low pressure region, the ACA simulation suggests 
that the navigation route should be shifted south of the Great 
Circle Route (Fig. 7). If the ship followed the Great Circle 
Route for which the distance is 4518.6 nautical miles, and 
travelled with the calm water speed of 18 knots, the navigation 


time under the weather conditions was 261.8 h (Fig. 8) and 
the average speed was 17.26 knots. Although the distance 
generated by the ACA (equal to 4571.6 NM) is longer than 
that of the Great Circle Route, the navigation time was only 
255.8 h and the average speed was 17.87 knots. Then 6 h and 
the corresponding amount of fuel can be saved compared with 
the Great Circle Route. For the west bound navigation, i.e. on 
the southern route, a ship encountered a head wind and head 
seas due to low pressure to the south of the Aleutian Islands. 
This casted an adverse effect on the navigation. The ACA 
algorithm therefore modified the route and caused it to deviate 
to the north off the Great Circle Route. The route entered the 
Bering Sea to the north of Aleutian Islands, and met favourable 
following seas and good weather, as shown in Fig. 9 and 10. In 
this case, the navigation time following the Great Circle Route 
was 267.2 h and the average speed was16.91 knots. The route 
produced by the ACA had a longer distance of 4571.72 NM, 
but the navigation time was only 264.7 h (refer to Fig. 11) 
and the average speed was 17.27 knots. The saving of 2.5 h 
and the corresponding amount of fuel was achieved. The 
simulation results shows that the weather routing by means of 
the ACA makes use of favourable meteorological conditions 
and avoids to some extent encountering head seas and heavy 
wind regions, reducing this way the risk of ship damage and 
cargo loss. As regards the execution time, although algorithm 
efficiency for the entire program has not been optimized, the 
average execution time was about 5 min, which is sufficient 
for decision-making in real-time navigation. 
261 
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Fig. 8. Iteration history of optimum navigation time for east bound 
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Fig. 6. Route simulation results for east bound navigation 
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Fig. 9. Route simulation results for west bound navigation 
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Fig. 10. Weather conditions and the simulated route for west bound navigation 
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6. CONCLUSIONS 


With the substantial rise in international oil prices, and global 
warming on the rise, to reduce operational fuel consumption 
and decrease air pollution has become one of the key goals of 
green shipping. In order to reach the goal, ships in operation 
can be fitted with energy-conserving equipment. However, 
optimum ship routing is the most direct and economical 
method. Therefore, application of decision support tools to 
help the routing of transoceanic crossings, not only can enable 
route planning automation, achieve the aim of carbon oxide 
reduction and energy conservation but also provide a reference 
for decision making in route planning. In this study a proper 
grid system and the improved ACA was used to establish 
a route searching model that simulates a living organism’s 
optimum behaviour. In addition, GA’s crossover and mutation 
computations, a ship’s critical speed function, and GIS’s 
spatial analysis and computation capability were integrated to 
automate the generation of a transoceanic crossing route, in 
consequence increasing route planning speed and quality, and 
lowering the ship crew’s burden. With the growing popularity 
of distributed computing environments, many applications of 
parallel processing and distributed processing have been widely 
implemented. Considering the efficiencies brought about by 
distributed processes, the using of the distributed parallel 
processing structure together with the feedback mechanism 
of the ACA collective intelligence can further enhance ACA 
resolution efficiency. And, the computational speed will be 
higher than that of other resolution methods. Therefore, in the 
future the size of the grid system can be reduced to conduct 
more refined simulations by means of the combination of more 
waypoints. Also, in this study the candidate routes that pass 
through obstacles was removed without integrating an obstacle 
avoidance mechanism (such as island avoidance navigation). 
A relevant coastal navigation routing technique, suggested by 
Tsou [23], could be added in the future so that the results might 
better satisfy requirements of real time navigation decision 
support and enhance the energy conservation effect. 
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ABSTRACT 


Digital signal processing technology has revolutionized a way of processing, visualisation and interpretation 
of data acquired by underwater systems. Through many years side scan sonars were one of the most widely 
used imaging systems in the underwater environment. Although they are relatively cheap and easy to deploy, 
more powerful sensors like multibeam echo sounders and sonars are widely used today and deliver 3° 
bathymetry of sea bottom terrain. Side scan sonar outputs data usually in a form of grey level 2° acoustic 
images but the analysis of such pictures performed by human eye allows creating semi-spatial impressions 
of seafloor relief and morphology. Hence the idea of post-processing the side scan sonar data in a manner 
similar to human eye to obtain 3° visualisation. 
In recently developing computer vision systems the shape from shading approach is well recognized 
technique. Applying it to side scan sonar data is challenging idea used by several authors. In the paper, 
some further extensions are presented. They rely on processing the backscattering information of each 
footprint (pixel in sonar image) along with its surroundings. Additionally, a current altitude is estimated 
from the size of shadow areas. Both techniques allow constructing 3° representation of sea bottom relief 
or other investigated underwater objects. 


Key words: side scan sonar data processing; shape from shading 


INTRODUCTION 


Acoustic sensors offer robust insight into underwater 
environments, as they are characterised by longer range and 
wider angle coverage compared to video cameras or other 
sensors. Also, they are capable of providing satisfactory 
results of mapping the environment in turbid waters. Side scan 
sonar (SSS) is one of the most widely used imaging systems 
in underwater environment. It is relatively cheap and easy to 
deploy, in comparison to more powerful sensors like multibeam 
sonars. However, it has some limitations, such as its inability 
to recover the seafloor depth or submerged object information 
directly [1]. 

Conventional side scan sonars are used until today. They 
emit fan-shaped pulses down toward the seafloor across a wide 
angle perpendicular to the path of the sensor through the water, 
which may be towed from a surface vessel or mounted on the 
ship’s hull. The intensity of the acoustic reflections from the 
seafloor of a fan-shaped beam is recorded in a series of cross- 
track slices. When stitched together along the direction of 
motion, these slices form an image of the sea bottom within 
the swath of the beam. 

During years of evolution side scan sonar were improved 
including higher frequency systems and the dual-frequency 


sonars, and the combined side-scan and sub-bottom profiling 
sonar. Generally, increasing frequency improves resolution of 
obtained images. The development of digital signal processing 
technology made possible real time processing of acquired 
data to form synthetic aperture. The principle of such synthetic 
aperture sonars (SAS) is to move a sonar along a line and 
illuminate the same spot on the seafloor with several pings. 
This produces a synthetic array equal to the distance traveled. 
By coherent reorganization of the data from all the pings, 
a synthetic aperture image is produced with improved along- 
track resolution. In contrast to conventional side scan sonar, 
SAS processing provides range-independent along-track 
resolution. At maximum range the resolution can be magnitudes 
better than that of side scan sonars. 

Until the development of so-called 3° side scan sonars the 
side scan sonar imagery were presented in the form of classical 
2D images. 3P side scan sonars uses interferometric acoustic 
technology and by application of digital signal processing 
techniques it operates over a field of view larger than 90° with 
no blind spots under the survey track. 

However, many SSS 2° images exist, that could be 
transformed into 3° representation in an algorithmic way using 
the echo intensity information contained in grayscale images. 
In the imaging sonar systems, the characteristics of acoustic 


POLISH MARITIME RESEARCH, No 3/2013 39 


Fig. 1. Sample sea bottom images obtained by three generations of side scan acquisition systems 
a) conventional side scan sonar b) synthetic aperture sonar c) 3° side scan sonar 


energy backscattering by a target localised on seafloor is 
utilised. Instead of arbitrarily using the formula like Lambert’s 
law to describe the seabed backscattering dependence on an 
incident angle [2], the authors propose their own approach to 
derive this dependence locally from experimental data. 

In this paper, two methods for 3° seafloor and submerged 
objects shape reconstruction from SSS images, based on 
shape from shading (SFS) approach, are presented. In the first 
method, for estimation of bottom local depth at a given pixel 
of sonar image, the information from both currently processed 
and previous ping is utilised. It allows the local seabed surface 
element orientation to have two degrees of freedom. Prior 
to the investigation of this approach on sidescan sonar field 
data, its verification on multibeam sonar backscatter data is 
presented, with use of the seabed 3° model (i.e. the set of (x, 
y, Z) seabed surface points derived from multibeam soundings) 
as a reference. 

The second method is used for 3P submerged object 
shape reconstruction and its imaging. The direct application 
of classical SFS technique seems not to be best suited to 
submerged object visualisation as it leads to obtaining very 
smooth shapes which differ significantly from actual forms of 
artificial objects. The proposed method, besides the utilisation of 
the local altitude gradient estimation by SFS algorithm with use 
of the backscattering coefficient angular dependence function 
derived experimentally, also estimates the elevation change on 
the basis of the dimension of acoustic shadow areas. 

Several techniques of 3P geometry reconstruction for 
seabed surface or submerged objects using side-scan sonar 
images has been reported [3, 4, 5]. Mainly, they use the 
techniques based on the problem inverse to image formation, 
namely Shape from Shading (SFS), which is one of classical 
problems in computer vision (see [6] for a collection of 
significant papers on SFS). The goal is to derive a 3P scene 
description from one or more 2° images. Given a grey level 
image, the aim is to recover the light source and the surface 
shape at each pixel of the image. 

Brief description of the SFS principles and basic methods 
along with literature references is given in [7]. 


SHAPE FROM SHADING FOR SIDE SCAN 
SONAR IMAGERY 


Seafloor relief reconstruction from sonar image 
In the construction ofa seabed elevation map from side scan 
images, the SFS technique relays on calculating the local slope 


of bottom relief, given the image pixel intensity, the assumed 
dependence of bottom surface backscattering coefficient 
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on incident angle (what corresponds to reflectance map in 
classical SFS), and the estimated local incident angle value. 
The presented algorithm for 3 seafloor relief reconstruction 
[8] has been developed assuming that the acoustic reflectivity 
model as well as the altitude H of the sonar transducer are given. 
The intensity (grey level) of a pixel in sonar image is treated as 
proportional to the acoustical intensity of backscattered echo. 
Also, the simple, straight line propagation path of acoustic wave 
in water column is assumed, and, the dimensions along vertical 
axis (z) of an object to be reconstructed are assumend to be 
small in comparison with the sonar transducer altitude. 

The geometry used in derivation of the reconstruction 
algorithm is presented in Fig. 2. The beam of a side scan sonar 
covers an angular sector from Q,,,, tO Omax 

The relation between the time instant t, in an echo envelope 
and the across-track coordinate x, of a corresponding point P, 
on seabed surface, may be expressed by: 


(1) 


where: 
t 2 2H/c 
c — sound speed in water. 

At the time instant t,, the seafloor surface S; is insonified, 
the area of which may be for a flat bottom case expressed by 
the classical equation [2]: 


CT 
UV isin 0; 2) 
where: 
©,- the along tract transducer beamwidth, 
R;— the range from the transducer to the point P.,, 
1 — the transmitted pulse length. 


The 3P bottom relief was reconstructed by estimation of 
an altitude z(x, y) sequentially for consecutive discrete points 
(x, y) on a plane, using the scheme depicted in Fig. 2. For 
the (i, j) iteration (where i — number of processed line in the 
sonar image corresponding to one sonar ping, j — number of 
pixel belonging to this line), i.e. the point P; = (x, Y; z(%,, yi) 
altitude estimation, the local triangle facet was being taken 
into account, with vertices at two previously estimated points 
Pi; = (% Yi Z(%, y,.1)) and P; ji = (x, Yo Z- yi), and 
currently estimated point P,. Using the applied model, the value 
chosen for z; allows for calculation of normal N;; to the surface 
facet, the angle ©, and the local intensity I, value, which then 
many be compared with that from the original sonar image. 


Ray to point 


Ray to point P} — Target plane 
Sonar kę = at point Pj 
motion haw Ni \ 
Sonar 


\ Seafloor surface Sij 
insonified at a given time tj 


PiE, yi, Z(Xj, y;)) 


Fig. 2. The geometry used in derivation of the seafloor relief reconstruction algorithm and its grid approximation 


Optimal z;, i.e. that giving I equal to a measured value, may 
then be estimated [7]. 

The reliable reflectance model (namely, here: the seafloor 
acoustic backscattering coefficient dependence on an 
incident angle) is a crucial element of the SFS-like 3P object 
reconstruction procedure. Some authors suggest that the 
Lambert’s Law provides a good fit to seabed backscattering 
dependence on an incident angle [2], but in general, the 
description of seafloor backscattering phenomenon is 
complicated. The backscatter from the seafloor is generally 
considered to be composed of a combination of surface and 
volume scattering, i.e. roughness interface scattering and 
scattering from inhomogeneities within the sediment volume. 
The detailed investigation of seafloor backscattering angular 
dependence, containing the theoretical modelling along with 
experimental verification, is described in [9]. It shows that 
although in several cases the application of Lambert Law might 
be allowed, in many others this law is not satisfied. This fact 
was also proved by further experiments. Taking it into account, 
the authors proposed their own approach here, which assumes 
the estimation of backscattering angular dependence coefficient 
locally for currently used experimental data. Assuming that 
we may choose the region within the processed sonar image, 
where the seabed surface is flat and composed of the material 
characterised by constant properties, we may state that the 
grey level of a given pixel in an analysed area should depend 
only on the experiment geometry (e.g. position of the source, 
the incident angle), the sonar calibration data, and the seafloor 
backscattering coefficient. If the investigated area covers some 
range of incident angle values, the backscattering coefficient 
angular dependence may be estimated for this range using pixel 
grey level values, as it is presented in the next section. 


The use of shadow area dimension for object 
height estimation 


One extension has been made for the proposed 3° shape 
reconstruction algorithm from SSS data when used to artificial 
objects lying on seabed (like wrecks for instance). Such object 
usually have regular shapes and their sonar images often are full 
of shadow areas. The proposed algorithm extension assumes 


the utilisation of the information which may be derived from 
shadow areas. Namely, the height of an object element is being 
estimated from the size of a current shadow area. Shadow area 
include the pixels from SSS image witch values are below 
chosen threshold intensity value I,,. In a case of a shadow 
zone detection of a length j along x axis, the altitude values 
from z,,, to Z; were set to unknown, and the z,,,,, value was 
calculated as: 


AX h 
tan Q; 


(3) 


Fig. 3 presents in a schematic way the influence of the I, 
on the algorithm results for two cases of 1„ values. Figure 3a 
shows the sample dependence of the intensity value on x 
coordinate for a given fragment of one line in side scan sonar 
image, along with two I,, values indicated. Figure 3b presents 
the reconstructed seabed altitude for these two cases. For l, 
case, the larger shadow zone occurs and no reconstruction is 
obtained inside [x,, x,] range. 


Zięj+1 7 Zi 
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Xa 
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Fig. 3. The influence of the shadow threshold value I, on the altitude 
reconstruction algorithm results 
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Fig. 4. Seabed backscattering coefficient angular dependence averaged for 40 multibeam soundings for each of 4 sites 


Fig. 5. Verification results of 3° shape reconstruction method using multibeam sonar data for 4 sites in Gdansk Bay region. 
First column: sonar images constructed using backscattering strength of echoes corresponding to several multibeam soundings 
(horizontal lines in an image) and beams (vertical lines in an image). Central column: visualisations of seabed 3° models derived 
from multibeam soundings by the EM3002 sonar software in the form of sets of (x, y, z) seabed surface points — these seabed surfaces 
were used as a reference data. Right column: 3° seabed surface shapes reconstructed using the proposed SFS-like algorithm 
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EXPERIMENTS AND RESULTS 


Experimental derivation of backscattering 
coefficient angular dependence 


As it has been introduced in the previous section, the seafloor 
backscattering coefficient angular dependence function, needed 
in the SFS-like 3P shape reconstruction procedure, has been 
derived experimentally by choosing flat surface area for this 
purpose. 

The data used in this experiment were acquired by the 
Kongsberg EM 3002 multibeam sonar in Gdansk Bay region 
of the Baltic Sea for 4 separate sites characterised by different 
seafloor type. The sonar operating frequency was 300 kHz, the 
width of beams: 1.5° x 1.5°, the transmitted pulse length: 0.15 
ms, the echo sampling rate: 14.3 kHz. The bottom depth was 
in a range approximately between 10 m and 100 m. 500 swaths 
from each of four seafloor types were processed. For each 
swath, 160 beams covered the angle sector from -65° to 65°. 

Fig. 4 presents the seafloor backscattering coefficient 
angular dependence derived by the described approach and 
averaged for 40 multibeam soundings for each of 4 sites. The 
significant differences in the results between several cases are 
well visible. It proves the expectation that this dependence 
characterises the seabed locally and therefore it should be 
preferable to derive it experimentally in-situ. 


Verification of the shape reconstruction 
procedure using multibeam sonar data 


The mutibeam data described in the precious subsection 
were also used in the verification procedure for the proposed 
algorithm of 3° shape reconstruction from acoustic data. 
The obtained results are presented in Fig. 5 for 4 different 
acquisition sites in Gdansk Bay. The seabed “sidescan sonar- 
like” images constructed from multibeam backscatter data 
(first column of pictures in Fig. 5) were processed by the 
proposed SFS-like algorithm to reconstruct the seafloor relief. 
The results are shown in the right column of pictures in Fig. 5 
while the reference data, i.e. 3° seabed surfaces derived from 
multibeam data by EM3002 sonar software, are shown in the 
central column in this figure. 

It is actually visible in Fig. 5 that although the proposed 
method introduces some “artifacts” in across-track direction 
(from bottom left to top-right on the pictures in Fig. 5), the 
differences between the reconstructed and the reference seabed 
surfaces are small enough to be primarily acceptable. For the 
quantitative estimation of the proposed method performance, 
the rms difference between the reconstructed and the reference 
surface has been calculated for each of 4 cases. The obtained 
rms values were between 3 cm for the first site (a) and 9 cm 
for the last site (d) what is a good result taking into account 
that the geometric dimensions of the experiment were of tens 
of meters. 


Application of the 3? shape reconstruction 
procedure to side scan sonar data 


The developed procedure of 3° seafloor relief reconstruction 
was also tested on DF1000 dual frequency (100 and 500 kHz) 
side scan sonar data collected in the Gulf of Gdansk. For 
further processing the 500 kHz data were used as they represent 
the higher resolution entity. The specific areas were selected 
for calculation of seabed backscattering coefficient angular 
dependence. The results of 32 seabed surface reconstruction 
are presented in Fig. 6. 


< y =200 m > 


«x=45m-> 


Fig. 6. The sample seabed image with selected areas used for seabed 
backscattering calculation [10] 


The developed procedure of object 3° shape reconstruction 
was tested on side scan sonar data downloadable from Marine 
Sonic Technology, Ltd. website. The sample car tire image 
acquired by side scan sonar is presented in Fig. 7 and its 
reconstruction is shown in Fig. 8. 


Acoustics 
shadow 
areas 


2m 


1—2 m—\—— 


Fig. 7. Side scan sonar image of car tire from sonar operated 
at 2400 kHz [10] 


Fig. 8. The 3° car tire reconstruction 


The next test was performed on side scan sonar data 
containing the image of USS Utah as shown in Fig. 9. The 
image region within the rectangle was used to verify the 
performance of the algorithm. 
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A 270 m 
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Fig. 9. Image of the USS Utah, resting in Pearl Harbor near Ford Island, acquired by the U.S. Army 
7 Engineer Detachment using Sea Scan Centurion system operated at 600 kHz [10] 


The comparison of the results obtained when applying two 
different shadow threshold values In, = 0.051 nax and Li = O. rax 
is presented in Fig. 10. It may be seen that for lower I,, value 
case (Fig. 10a), more details of the reconstructed shape may 
be visible than for higher I, value (Fig. 10b). But on the other 
hand, if the I, is too low, the artefacts may occur due to using 
the information from very dark, and possibly noised pixels, for 
the object local altitude slope estimation. 


Fig. 10. 32 wreck shape reconstruction results using two different values of 
shadow threshold 


CONCLUSIONS 


The methodology based on shape from shading approach 
for 3° reconstruction of seafloor and submerged objects from 
side scan sonar records was presented. It was verified using the 
data acquired by multibeam sonar system that contains both 
backscatter and bathymetry data and also tested on side scan 
sonar data. The principle advantage of the methods is their 
simplicity and the ability to produce the results within sequential, 
one-run processing. It was also shown that angular dependence 
of backscattering plays a significant role during image processing 
and should be estimated for a observed subareas. 

At the current stage of the development the applicability 
of the presented solution might be questionable and further 
investigation are required. In particular, the authors expect 
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that presented seafloor reconstruction algorithm performance 
can be improved by applying procedure similar to synthetic 
aperture sonar principle. The implementation of more 
advanced shape from shading algorithms could be required 
in more complex cases. 
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ABSTRACT 


The main goal of this paper is to describe the results of sounding the Gdansk Bay sea bed by using a SES- 
2000 Standard parametric sub-bottom profiler. Quality of data obtained during trials depends inter alia 
on a proper location of transducer array to reduce influence of pitch, roll and heave motions as well as 
ship self noise (resulting from bubbles due to propeller and flow around hull, vibration generated by main 
engine and auxiliary devices). Furthermore, calibration of complementary units such as GPS, heading 
sensor, MRU-Z motion sensor and navigation devices make sea-bed investigating system capable of working 
with its whole capability. Results of so prepared surveys have been presented and discussed. They contain 
also an elaborated map of Gdansk Bay with preliminarily classified sea-bed materials and description of 
most interesting areas. 


Key words: sea bed sounding; parametric sub-bottom profiler; sea bed structure 


INTRODUCTION 


Examination of the upper layer of sea bed can be performed 
with the use of invasive methods as well as non-invasive ones. 
The invasive methods apply only to the local structure of the 
sea bed, and are associated with geographical coordinates of the 
bore-hole point (point of evaluation). In addition, in such case 
a subtle structure of actual image of the bottom is destroyed. 
The non-invasive methods are mainly seismic ones which 
are rather suitable for determining the geological structure of 
the sea bed in relation to large depth [14, 16]. However, for 
the bottom layer they are relatively inaccurate. An alternative 
solution is to use acoustic methods which make accurate 
determining the sea-bed stratification, possible. These methods 
are non-invasive and relatively versatile. However there is 
a problem of penetration of bottom sediment by acoustic waves 
[1, 2, 9]. Attenuation of elastic waves by bottom sediments is 
about 1000 times larger than by the water [13]. For this reason 
the high-frequency waves are mostly reflected upon the water 
- bottom border surface. 

With taking into account the relevant constraints a parametric 
source of acoustic waves can be applied in such case. The 
principle of operation of such source is that, as a result of the 
nonlinear interaction between two large-intensity acoustic 
beams in the area of common impact of the beams, a new 
wave is formed in the area of a volumetric character [12]. As 
the beams called primary, differ only slightly in frequency, 
we get a wave called the differential wave whose frequency is 


equal to the difference in frequency of the primary waves. Such 
wave maintain the features of primary waves and therefore is 
of a high directivity. Due to the area of its generation which 
is three-dimensional, it has no additionally adverse side lobe. 
The properties were used in the construction of the parametric 
sonar, inter alia in the SES 2000 Standard parametric echo- 
sounder [17]. Extensive research of sea bottom in Gdansk 
Bay waters was performed with the use of this echo sounders 
[3-8, 10, 11]. 


RESEARCH VESSEL 


The platform for sea trials was set on board s/y Windspeel, 
a small research vessel. The main device was the SES-2000 
Standard sub-bottom echo sounder which consisted of the main 
unit situated below the deck and the transducer array mounted 
on aluminium arm starboard. The unit makes it possible to 
radiate low frequency sounding pulses set by the user within 
the range from 4 to 15 kHz and a definite changeable length. An 
interesting function during sounding was the option of obtaining 
sub-bottom data by applying two or three different pulses with 
user - defined frequencies. During the sea trials in question the 
frequencies of 4 kHz, 8 kHz and 12 kHz were set and combined 
to each other. As the best weather conditions are during calm 
sea which is not always achievable, hence a beam steering to 
compensate for vessel pitching and rolling was necessary. 

The obtained data were stored on hard disc in the two 
formats: *.SES one which includes information about signal 
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envelope and *.RAW data containing full waveform. The whole 
unit was fed with 230V (AC) provided from DC-AC converter 
connected to the set of three 110 Ah gel batteries charged 
during trials from solar batteries and a generator driven by 
vessel’s engine. During breaks in measurements the charging 
was provided by on-board 230 V Honda electric generating 
set. Additionally sea bed was sounded with the use of the 
multi-beam EM3002 echo sounder whose transmitted pulses 
were synchronized by means of an external trigger installed 
in the SES-2000 sounder that ensured undisturbed results. 
Complementary devices like the MRU-Z sensor of pitch, roll 
and heave motions precise GPS with heading sensor were 
installed and calibrated on the vessel. Signals from the sensors 
were distributed to the devicesby means of RS-232 splitters 
During the trials sound velocity profiles were stored in SD204 
unit and put in to controlling programs of sub-bottom and 
multi-beam echo sounder. The investigated area and transects 
selection along which the sounding took place were determined 
on the basis of geological maps delivered by Panstwowy 
Instytut Geologiczny (State Geological Institute) [15]. The 
appropriate data were delivered to the navigation software 
(Nobeltec) connected to autopilot, that ensured a satisfying 
accuracy of sounding along the set tracks. To help the skipper 
to control situation on sea two monitors were duplicated in 
cockpit to display navigation map and sounding results. The 
research vessel during measurement trials on Gdansk Bay is 
presented in Fig. 1. 

The devices necessary to ensure proper work of the whole 
unit were placed on specially prepared working racks inside the 
vessel (Fig. 2). Before beginning the trails on sea the units were 
calibrated in harbour and checked during work with different 
parameters and configurations. Additionally, the electric system 
was tested and checked against possible emission of noise 
which could influence sounding equipment. 

Special RS232 splitters were implemented into the system 
to avoid duplication of the units delivering information such 
as position and motion data. 


OBTAINING SEDIMENT SAMPLES 
During all sea trials a diver prepared for such action was 


collecting samples of material taken from bottom’s surface and 
other layers according to a type of sediment. The samples were 
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Fig. 2. Complementary devices 


stored in special containers, and then examined and classified. 
Additionally, it was necessary to collect enough amount of 
materials to be used during trials in 30 m long water tank of 
Gdansk University of Technology, where influence of different 
kind of materials on the parametric echo sounder records was 
measured and then used to prepare database for comparison of 
results from sea trials with data obtained during the “artificial 
sounding”. 


Fig. 3. Diver who prepares for taking bottom sediment samples 


Fig. 4. Mud dug out of sea bottom 


An example container filled with mud taken from sea 
bottom is presented below. Collecting such sediment was most 
difficult for the diver because of no visibility due to cloud of 
silt produced during digging the material. 

Enlarged pictures of sediments are presented in Fig. 5 
and 6. The example clearly shows difference in bottom layer 
structure. The example IV illustrates the sample of sand which 
results in high power echo compared with the sample ofmud 
of a low penetration - example III. 

Owing to the data which uncover clear parts of the cliff 
it was possible to elaborate pictures of the structure of such 
sediment layers by collecting the sediments and comparing 
them with the samples taken from sea bottom. Real image of 
such structure provides more information on its geological 
nature. In Fig. 7 is presented an example of enlarged pictures 
of materials taken from the pointed parts of the cliff material 
layers. 
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Fig. 6. Sediment samples in macro view, part 2 
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Fig. 8. Map of Gdansk Bay which shows classification of the bottom 
penetration areas. A - the area called “hilly sub-bottom”’; B - the area 
of medium sounding depth (up to 5 m penetration); C - the area of the 
deepest penetration (up to 30 m); D - the area of a highly reflective bottom 
(no penetration); E - the area of a mixed type of penetration (from 1 m up to 
10 m); Thick red line - a part of the transect E-F - taken for comparison 


EXAMPLES OF SEDIMENT 
CLASSIFICATION 


After several trails in Gdansk Bay it was possible to elaborate 
a basic classification of sub-bottom penetration areas indicated on 
the map (Fig. 8) which shows five main types of the results. 

Sample echograms corresponding to the marked areas are 
presented below: the first (Fig. 9) called “hilly sub-bottom”, 
shows a typical bottom which reflects strong echo, and another 
bottom profile at first glance looking like “hills”, where the top 
of it sometimes touches the bottom and even penetrates it e. g. in 
the point I. After examination, sediments of sand silt and clayey 
silts were found, what confirmed information taken from the 
geological map. From the point of view of the layered structure 
the area E is rather less interesting as it combines places with 
penetration up to 10 m with those of strong echo only and of 
no information of extra layers which corresponds to sediments 
like sand from which whole sounding pulse is reflected. The 
looking and analyzing of such long echograms gives possibility 
to determine permanence of sub-bottom layer. 

The other results were also compared with the data delivered 
by Panstwowy Instytut Geologiczny, representing a geological 
map of bottom structure (Fig. 11), prepared by using different 
methods including taking deep sediment samples. When 
analyzing echogram results shown in Fig. 12 one can notice 
that a part of sounding which dealt with sand characterized 
shallow penetration, and in the places where silt and clay 


8 m penetration 


Fig. 10. The area of a mixed type of penetration depth (from 2 m up to 10 m) — area E 
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Fig. 11. Geological structure of the marked transect 
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Fig. 12. Echogram taken along the transect marked with red line in Fig. 11 


occurred deeper souding results were obtained. Comparison of 
the results was possible due to the use of devices to precisely 
contol of the research vessel motion and sail along selected 
transects. During the measurements the parametric echosounder 
was tested with different sounding pulse frequencies to acquire 
information about echosounder parameters during sounding the 
deep bottom penetration areas. 

The last results come from the mixed areas B, C and D. They 
were taken along the transect marked with red line depicted on 
the geological map. 

When analyzing of the results (Fig. 14) is started from the 
left side, the data obtained in the area B representing strong echo 
reflected from hard bottom are there visualized. When passing 
to B sector a deeper penetration can be noticed - not as deep 
as in the middle of this area but with the same characteristic 
structure. The end of the echoogram shows rapid change of 
penetration depth compared with bathymetric data and change 
È of colour of the surface from orange to red which is related to 
Fig. 13. Geological map with the marked transect a higher level of the bottom echo strength. 
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Fig. 14. Echogram taken along the transect marked on the geological map (Fig. 13) 


SUMMARY 


In the paper is presented an experimental acoustic method 
for determination of sea bed material as well as structure of 
the upper layer of the sea bottom. 

After performance of trials and data processing, their results 
were analyzed and compared with those obtained during the 
trials in the water tankthat was satisfactorily implemented into 
the first step of classification of sediment layers. 

Additionally experiments in preparing long echograms 
(over a few kilometers) clearly showing changes in sediment 
structure, have ended with good results that was confirmed by 
the promising spectrograms presented in this paper. 

The systematic approach to the problem, consisted in 
starting with measurements of some acoustic parameters of 
sediments in laboratory conditions and then passing to the 
in situ investigations, resulted in determination of different 
properties of sea bottom in some areas of Gdansk Bay. 
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The Ship Handling Research and Training Centre at Ilawa is owned by the Foundation for Safety of Navigation 
and Environment Protection, which is a joint venture between the Gdynia Maritime University, the Gdansk University 
of Technology and the City of Ilawa. 


'Two main fields of activity of the Foundation are: 


O Training on ship handling. Since 1980 more than 2500 ship masters and pilots from 35 countries were trained at 
Iława Centre. The Foundation for Safety of Navigation and Environment Protection, being non-profit organisation 
is reinvesting all spare funds in new facilities and each year to the existing facilities new models and new training 
areas were added. Existing training models each year are also modernised, that's why at present the Centre represents 
a modern facility perfectly capable to perform training on ship handling of shipmasters, pilots and tug masters. 


Research on ship's manoeuvrability. Many experimental and theoretical research programmes covering different 
problems of manoeuvrability (including human effect, harbour and waterway design) are successfully realised at 
the Centre. 


The Foundation possesses ISO 9001 quality certificate. 
Why training on ship handling? 


The safe handling of ships depends on many factors - on ship's manoeuvring characteristics, human factor (operator 
experience and skill, his behaviour in stressed situation, etc.), actual environmental conditions, and degree of water 
area restriction. 


Results of analysis of CRG (collisions, rammings and groundings) casualties show that in one third of all the 
human error is involved, and the same amount of CRG casualties is attributed to the poor controllability of ships. 
Training on ship handling is largely recommended by IMO as one of the most effective method for improving the 
safety at sea. The goal of the above training is to gain theoretical and practical knowledge on ship handling in a wide 
number of different situations met in practice at sea. 


For further information please contact: 
The Foundation for Safety of Navigation and Environment Protection 


Head office: Ship Handling Centre: 
36, Chrzanowskiego Street 14-200 IŁAWA-KAMIONKA, POLAND 
80-278 GDAŃSK, POLAND tel./fax: +48 (0) 89 648 74 90 
tel./fax: +48 (0) 58 341 59 19 e-mail: office@ilawashiphandling.com.pl 

e-mail: office@portilawa.com 
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GDANSK UNIVERSITY OF TECHNOLOGY 


is the oldest and largest scientific and technological academic institution in 
the Pomeranian region. The history of Gdansk University of Technology is 
marked by two basic dates, namely: October 6, 1904 and May 24, 1945. 


The first date is connected with the beginning ofthe technical education at 
academic level in Gdansk. The second date is connected with establishing of 
Gdansk University of Technology, Polish state academic university. Gdansk 
University of Technology employ 2,500 staff, 1,200 whom are academics. 
The number of students approximates 20,000, most of them studying 
full-time. Their career choices vary from Architecture to Business and 
Management, from Mathematics and Computer Science to Biotechnology 
and Environmental Engineering, from Applied Chemistry to Geodesics and 
Transport, from Ocean Engineering to Mechanical Engineering and Ship 
Technology, from Civil Engineering to Telecommunication, Electrical and 
Control Engineering. Their life goals, however, are much the same - to 
meet the challenge of the changing world. The educational opportunities 
offered by our faculties are much wider than those of other Polish Technical 
universities, and the scientific research areas include all of 21st Century 
technology. We are one of the best schools in Poland and one of the best 
known schools in Europe — one that educates specialists excelling in the 
programming technology and computer methods used in solving complicated 

scientific, engineering, organizational and economic problems. 


THE FACULTY OF OCEAN ENGINEERING 
AND SHIP TECHNOLOGY 


The Faculty of Ocean Engineering and Ship Technology (FOEST) as 
the only faculty in Poland since the beginning of 1945 has continuously 
been educating engineers and doctors in the field of Naval Architecture 

and Marine Technology. 


The educational and training activities of FOEST are supported by 
cooperation with Polish and foreign universities, membership in different 
international organizations and associations, as well as participation in 
scientific conferences and symposia. Hosting young scientists and students 

from different countries is also a usual practice in FOEST. 


The activities of Faculty departments are related to: mechanics and 
strength of structures, hydromechanics, manufacturing, materials and system 
quality, power plants, equipment and systems of automatic control, mostly 

in shipbuilding, marine engineering and energetic systems. 


FOEST is a member of such organizations like WEGEMT; The 
Association of Polish Maritime Industries and the co-operation between 
Nordic Maritime Universities and Det Norske Veritas. The intensive teaching 
is complemented and supported by extensive research activities, the core 
of which is performed in close collaboration between FOEST staff and 
industry. We take great care to ensure that the applied research meet both 
the long term and short term needs of Polish maritime industry. FOEST 
collaborates with almost all Polish shipyards. Close links are maintained 
with other research organizations and research institutions supporting the 
Polish maritime industry, such as Ship Design and Research Centre and 
Polish Register of Shipping, where several members of the Faculty are also 

members of the Technical Board. 


The Faculty of Ocean Engineering and Ship Technology is a unique 
academic structure, which possesses numerous highly qualified and 
experienced staff in all above mentioned specific research areas. Moreover, 
the staff is used to effective co-operation and exchange of ideas between 
specialists of different detailed areas. This enables a more integrated and 
comprehensive treatment of research and practical problems encountered 
in such a complicated field of activity as naval architecture, shipbuilding 

and marine engineering. 


The staff of the Faculty has strong international links worldwide, being 
members or cooperating with international organizations like International 
Maritime Organization IMO, International Towing Tank Conference ITTC, 
International Ship and Offshore Structures Congress ISSC, International 
Conference on Practical Design of Ship and other floating Structures PRADS 

just to name a few. 
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e-mail: sekoce(dpg.gda.pl 


52 POLISH MARITIME RESEARCH, No 3/2013 


Gdansk University 
of Technology 


Faculty of, Ocean woot 


Engineering w 
andShip JEGHNOJOGY. 


ZN ——— yh 


c 7 


www.oce.pg.gda.pl 


